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I. INTRODUCTION 


For many years the properties of hard x-ray sources in solar flares have been 
studied in order to gain greater insight into the underlying processes causing these 
highly energetic events. To date, temporal, spectral, and polarization measurements 
have failed to provide unambiguous information on the processes responsible for 
accelerating protons and electrons. There are a variety of reasons for this, foremost 
among which is the limited temporal and/or spatial resolution available. (Tandberg- 
Hans sen and Emslie, 1988). 

High resolution spatial observations show great potential for advancing our 
knowledge of the hard x-ray production mechanism since for the first time our 
ability to measure spatial structure approaches temporal and spatial resolution 
scales of physical interest (e.g., the collisional stopping length of a high-energy 
electron in the corona; Campbell, Davis, and Emslie, 1991). 

Unfortunately, hard x-rays (e > 10 keV) cannot be imaged by conventional 
methods such as lenses or mirrors. In addition, the usefulness of grazing incidence 
telescopes in this energy regime is limited due to the stringent requirements im- 
posed on the smoothness of the optics, the physical sizes necessary for a reasonable 
collecting area, and the unreasonably long focal lengths which would be required. 
On the other hand, images can be obtained in a number of ways from absorptive 
grids or masks made from high Z materials such as tungsten. 

One popular approach in the astrophysics community has been coded aperture 
telescopes (Murphy, 1990). These instruments use a multiple pinhole mask to 
produce a pattern which is the correlation function of the source image with the 
pinhole array. The source image can be recovered by inverting the correlation 
integral. Unfortunately, while this instrument allows 50% quantum throughput, 
practical limitations to its spatial resolution (i.e., = 25 arc seconds based in part 
on attainable detector pixel sizes) make it marginal for use in imaging flares. 


Another type of instrument, the Fourier telescope (see Chapter IV), offers in 
principle much better resolution but the throughput will be reduced to 25%. These 
instruments work by measuring the intensity distribution in Fourier space and then 
either (a) numerically recomposing the components or (b) “back projecting” the 
detector response onto the sky to provide an image (Murphy, 1990). 

Essentially the concept involves sampling selected two-dimensional Fourier 
components from a wavefront emitted by the source. Each component may be 
visualized as a single point on a common complex surface. By measuring a sur- 
prisingly small number of components over a sufficiently large spatial frequency 
spectrum, this Fourier surface may be approximated. An inverse Fourier transform 
of this surface function then yields a noisy or a dirty image (using terminology taken 
from radio astronomy). This description is aesthetically correct, for when viewing 
a typical dirty image, the processing artifacts can make the picture very difficult for 
the eye to understand. Techniques taken from radio astronomy such as the max- 
imum entropy and the clean methods can be used to suppress these artifacts and 
produce a more meaningful image (Campbell et al, 1991a). This will be discussed 
in greater detail in Chapter IV. 

Several Fourier telescope concepts have been developed to study solar flares. 
They include the High Energy Imaging Device (HEIDI), an imaging experiment for 
the Max ’91 balloon program planned for a 1992 launch (Figure 1.1); the Japanese 
Solar- A Hard X-ray Telescope on the YOHKOH satellite launched in August 1991; 
the High Energy and Spectroscopy (HEISPEC) investigation which is proposed as 
the key component of the High Energy Solar Physics (HESP) mission (Figure 1.2) 
scheduled for the next solar maximum; and the lunar Pinhole/Occulter Facility 
(P/OF). These instruments are expected to be able to provide full Sim coverage 
within the instrument field of view with resolutions approaching 1 arc second. 

Figure 1.3 illustrates the lunar P/OF. This is an innovative idea which takes 
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Figure 1.2: The HESP Concept 


A conceptual sketch of HESP, a hard x-ray imaging system mounted upon a free 
flying satellite platform. In this case the grid planes would be separated by 5 meters. 


4 


advantage of the relatively weak gravity of the moon to construct a large extended 
telescope structure for obtaining high resolution hard x-ray images. The larger ring 
rotates in the horizontal plane a full 360° while the small upper grid tray is free 
to move in a constrained path along the larger ring through a 180° extent. The 
pedestal at the center also supports a lower grid tray which is moved independently 
so as to maintain a parallel orientation and alignment with the upper one. This 
allows the observatory to track a point in the lunar sky (e.g., a solar flare) to within 
the accuracy of the pointing system (Campbell et al., 1991b). 

Figure 1.4 illustrates another concept in which a Fourier telescope (or a redun- 
dant coded aperture) is deployed from the shuttle. In this concept, the successful 
heritage of actual shuttle experiments employing extended truss structures in mi- 
crogravity would be exploited. During launch, the truss structure would be stored 
in a cannister within the shuttle bay. Once on orbit, the truss would be extended 
to its full length (e.g., 100 feet) and stabilized using small, microthrust nozzles in 
conjunction with a laser feedback system for control stability and pointing accuracy. 
For a Fourier telescope, one grid plane would be mounted at the top of the boom 
with the lower grid plane, detector arrays, and electronics mounted in the bay. 

Figure 1.5 illustrates how a free flying Fourier telescope could be employed us- 
ing extendable boom technology. Another approach would be to use an extendable 
tube to separate grid planes. The free flyer is probably the most cost effective of 
the options considered here in that it holds the potential for providing maximum 
science return for several years for a relatively moderate price. 

A. Research Objectives 

The primary objective of this research is to derive the hard x-ray emission as 
a function of time, energy, and loop position for two prominent solar flare models 
found in the literature and here denoted the Model T and the Model N-T. Once 
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Figure 1.3: The LUNAR P/OF 

A conceptual sketch of the lunar P/OF, a hard x-ray imaging system as part of 
a lunar-based facility. This concept would take advantage of the reduced gravity 
with an innovative pointing and tracking design. 
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Figure 1.5: A Fourier Telescope For A Free Flying Satellite 
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the haxd x-ray emission profiles axe obtained, distinguishing signatures may be 
identified. This will be discussed in detail in Chapters II and III. 

The next step will be to construct two numerical models representing the two 
basic classes of Fourier telescopes, a Spatial Modulation Collimator (SMC) and a 
Rotating Modulation Collimator (RMC). Each will view the previously obtained 
profiles to determine over what energy range, and at what integration times, good 
images can be obtained. Chapters V and VI will discuss these telescopes in detail. 
In addition, various instrument- related questions will be addressed such as low flux 
level, effective field of view, twist, and noise capabilities and limitations. 

B. Flare Background 

The Sun’s atmosphere consists of four regions each defined by different char- 
acteristics. The lowest is called the photosphere and is mostly what we observe 
when we look at the Sun in visible light. It is (relative to the other layers) dense 
and opaque and on the order of 500 km in thickness. Temperatures in this region 

are on the order of 6000 K and densities are 10 1 cm 

Above the photosphere lies the chromosphere, which is relatively more trans- 
parent and about 2500 km thick. Temperatures are on the order of 5000 K and 
densities are typically 10 12 cm" 3 . Above the chromosphere lies the transition 
region which is on the order of 200 km thick in the quiet Sun. In a flare, this 
thickness may decrease to less than 1 km. In this region, the temperature rises 
rapidly approaching 10 6 K at the bottom of the next region which is called the 
corona. The corona may extend out for thousands of kilometers and is charac- 
terized by relatively high temperatures and low densities (i.e., < 10 9 cm" 3 ). In a 
flare, the density within the corona may increase to 10 11 cm -3 . 

Flares occur primarily in active regions on the Sun. These active regions 
extend vertically from deep photospheric levels through the chromosphere and the 
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transition region high into the corona. They exhibit much higher magnetic field 
strengths (i.e., 100-1500 G) than the surrounding quiet regions (1-10 G). Soft x-ray 
observations from Sky lab showed that to a large extent the active solar plasma is 
contained in loops that exist on scales from the smallest 1 arc second diameter 
bipolar loops emerging from below the photosphere to enormous coronal loops and 
arches spanning distances of a solar radius or more. 

A typical flare evolutionary sequence starts with precursors, such as gradual 
soft x-ray brightenings, which lead into an impulsive phase which, in turn, is fol- 
lowed by a gradual phase and then by a decay. It is during the impulsive phase that 
most of the energy which has slowly been built up and stored in the local magnetic 
field during the pre-flare period is suddenly released. Observations show this phase 
to be characterized by intense, rapidly fluctuating bursts of high-energy radiation, 
such as hard x-rays and 7— rays; bursts of accelerated electrons and protons; and 
plasma and bulk fluid acceleration. 

There are also associated emissions at other energies, such as EUV and opti- 
cs* generated as a result of the thermal response of the atmosphere to this phase. 
Figure 1.6 compares microwave and hard x-ray signatures for impulsive flares, es- 
tablishing the simultaneity between the two. Since the microwaves are produced 
by the gyrosynchrotron radiation from the hot electrons spiraling around magnetic 
field fines in the loop, this figure illustrates the importance of complementary mea- 
surements and the importance of microwaves as a diagnostic tool. Figure 1.7 shows 
the relationship between two emission fines and hard x-rays for a selected flare and 
illustrates the importance of complementary x-ray and UV observations. In partic- 
ular, UV observations can provide information about preheating processes prior to 
the impulsive hard x-ray bursts. These figures also are examples of the impulsive 
phase of a flare. 

In order to produce hard x-ray and 7-ray bursts, the basic flare mechanism 
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Figure 1.6: Typical Microwave And Hard X-ray Emissions 


These observations indicate that both the microwave and hard x-ray emissions 
from the impulsive phase of the flare occur simultaneously (Kosugi, Dennis, and 

Kai, 1988). 


11 




TIME (U.T.J 

Figure 1.7: Selected Emission Lines And Hard X-ray E mis sion 

This shows the behaviour of two ultraviolet (oxygen and iron) lines overlaid on the 
impulsive phase x-ray and microwave emissions. The UV lines should be a good 
impulsive phase diagnostic since it provides information about preheating processes 
in the flare prior to the onset of the impulsive phase. T his illustrates the importance 
of hard x-ray and UV complementary observations. Time profiles are from the 
impulsive phase of a flare observed on November 1, 1980 (from Tandberg-Hanssen 
and Emslie, 1984). 


12 



must be capable of accelerating large quantities of particles, both electrons and 
protons. This particle acceleration is also implied by centimeter microwave bursts 
and energetic proton streams. Microwave bursts, while relatively low in energy 
compared to the hard x-ray bursts, provide significant insight into the processes at 
work. Proton streams at the Earth typically have a power law shape (van Hollebeke 
et al., 1975). The accelerating mechanism is unknown and continues to receive a 
great deal of attention from the scientific community. 

Hard x-ray emission in flaxes can range from the simple isolated bursts with a 
total duration of several seconds to extended bursts lasting up to 10 3 sec, such as 
the great flare of August 4, 1972 (Hoyng et al., 1976). Figure 1.8 illustrates a hard 
x-ray burst time profile and indicates that the same mechanism produces photons 
over a the range from 28-111 keV. In addition, the rapid falloff of intensity as a 
function of energy is shown. Other observations have extended this idea to even 
higher wavelengths. This figure also indicates that this flare is composed of several 
events occurring somewhat sequentially in time. 

Simple bursts are characterized by a few loops sharing common footpoints 
while extended bursts are characterized by arcades of loops extending over the 
size scale of an active region (i.e., 1 arc minute). In both cases, strongly sheared, 
magnetic field lines running almost parallel to the neutral line, occur in this region. 
One scenario for energy release involves a long filament or strand of solar plasma 
several arc seconds long erupting into the tangle of stressed loops triggering an 
energy release in one or more. As individual loops are triggered they in turn trigger 
adjacent loops with some time delay. The flare then becomes a complicated cascade 
of sequentially triggered loops or elementary bursts. 

Figure 1.9 shows an M class flare observed by the author at the Marshall 
Space Flight Center Solar Magnetograph Observatory in H-alpha and illustrates 
the complexity of these events when observed in H-alpha. Clearly, complementary 
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Figure 1.8: Time Dependence Of A Typical Hard X-ray Burst 

This shows the simultaneity in three different energy bins of hard x-rays from 28 
to 111 keV indicating that the same mechanism produces the x-ray photons over 
these energies. In addition, the rapid decrease in intensity as a function of energy 
is shown, (from Kane et al., 1979b). 
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observations will be helpful in providing greater insight. 

Peterson and Winckler first detected hard x-rays from solar flares in 1959 
using a balloon observatory. Since that time, space and ground-based instruments 
have expanded the observational data base across the electromagnetic spectrum 
including hard x-rays presumably produced by bremsstrahlung. X-ray emission by 
bremsstrahlung (braking radiation ) is simply described as photon generation due 
to the large angle collision and the subsequent deceleration of an energetic electron 
with an ion in which the electron remains free. 

Understanding hard x-rays is extremely important to understanding the over- 
all physical processes occurring in the flare since they may be directly related to 
the dominant energy release processes occurring during the event. For example, 
bremsstrahlung resulting from electron beam deceleration in a thick target (see 
Chapter II) is extremely inefficient (= 10 5 erg of electron energy is needed to pro- 
duce = 1 erg of x-ray energy; Tandberg-Hanssen and Emshe, 1988). Consequently, 
extremely high limits are placed on the minimum required energy m the release. 

Not only is the energy required in such a model large, it must be in the low 
entropy form of directed, accelerated particles rather than a high entropy release 
in the form of bulk heat (Smith, 1980). Thus, an unambiguous verification of the 
nonthermal model would provide tremendous insight into the energy release process 
and its demands on the magnetic field. 

These constraints make attractive a model in which the bremsstrahlung pro- 
ducing electrons do not interact with cooler electrons thereby allowing a more effi- 
cient conversion of their energy into x-rays. Such a model could in principle reduce 
the number of electrons required at the beginning of the release (Tandberg-Hanssen 

and Emslie, 1988). 

Such a thermal model, the Model T, may be pictured simply as a hot mass 
of plasma with temperatures sufficiently large that 
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Figure 1.9: M Class Flare In H-Alpha 

This indicates the complexity of even moderately sized flares as seen in H-alpha. 
The flare was observed by the author on June 1, 1991 at the Marshall Space Flight 
Center Solar Magnetograph Observatory. 
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lcT « e « 10 fceV 


whence 


T « 10 8 K. 

However, one runs into difficulties with this picture as well. Kahler ( 1971a, b) 
pointed out that at the 10 8 K and upward temperatures required of such a thermal 
hard x-ray source, conductive cooling (using classical relationships, Spitzer, 1962) 
would be so efficient that a supply of energy even larger than that required for 
a nonthermal model would need to be provided for the observed duration of the 
burst. He further pointed out that it would be exceedingly difficult to confine such 
a hot plasma because of the very long collisional mean free paths 


A 


2 X 10 4 T 2 10 20 ^in9 
« > 10 

n n 


cm 


of the electrons and the consequent high probability of escape. 

Brown et al. (1979) showed, however, that collective plasma processes can 
significantly alter the physics of such sources. A very large electrical current would 
be generated as a result of the electrons streaming out of the source into the sur- 
rounding cooler plasma. The electric field set up by the resulting charge separation 
would accelerate ambient electrons into a reverse current with a velocity relative to 
the background ion distribution sufficiently large that ion-acoustic waves would be 
generated. These waves would be turbulent in nature and would tend to form into 
a barrier region which would reduce the mean free path of the electrons resulting 
in more effective confinement (Brown, Melrose, and Spicer, 1979). 

The turbulent containment region would quickly evolve to a marginally stable 
state (Manheimer, 1977) in which the escape of the most energetic electrons would 
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be balanced by a return current just sufficient to maintain the barrier. The region 
containing the turbulence would be relatively thin (i.e., on the order of 1 km) 
and would move along the loop at the ion sound speed driven by the temperature 
gradient. Chapter III will discuss this in more detail. 

Unfortunately, both of these models predict similar hard x-ray spectral signa- 
tures (Brown and Emslie, 1987; Emslie and Machado, 1987; Emslie, 1989). For 
example, Figure 1.10 shows an example of the hard x-ray spectra for four typical 
flares spanning a range from 10 to 100 keV (Kane et al., 1979b). These curves 
may be fitted by a power law, however this is not necessarily a unique mathemat- 
ical solution. The fact that these observations alone are insufficient to distinguish 
between models serves as a motivation to move toward high resolution, imaging 
telescopes which can resolve directly the evolution of hard x-ray sources in flares in 
order to distinguish between these scenarios. 

As shown in Figure 1.11, Model N-T can be qualitatively characterized by an 
initial footpoint structure, which subsequently evolves into a filled loop. Chapter II 
will discuss this in greater detail. By contrast, Model T initially exhibits a bright 
central region in addition to footpoint emission formed by high energy electron 
escaping from the central region. This central region grows and gradually fills the 
loop. Chapter III will discuss this in greater detail. These two different spatial 
emission signatures may provide the basis for a discriminating test between the two 
models through imaging (Emslie, 1981). 

However, this probably will only be true at the beginning of the imp uls ive 
phase. At the end of the impulsive phase (typically 30-60 sec), the loop should 
become uniformly filled in either case and if viewed at this time would probably 
offer little insight into the mechanism at work. In addition, since the impulsive 
phase itself is composed of many elemental flare bursts, it may prove necessary to 
image the first one as subsequent bursts may overlap and obscure signatures. 
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Figure 1.10: Hard X-ray Spectra For Four Different Events 

This shows four typical observations spanning a spectral range from 10 to 100 keV 
(Kane et al., 1979b). Note that each observation can be fitted by a power law; 
however, this may not be a unique solution. These observations alone axe not 
sufficient to provide an unambiguous determination of impulsive phase processes. 
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Direct imaging of hard x-ray sources in solar flares has only been accomplished 
fairly recently with the Hard X-Ray Imaging Spectrometer (HXIS) instrument on 
the Solar Maximum Mission (SMM) spacecraft and with the Solar X-Ray Telescope 
(SXT) onboard the Hinotori satellite. Hoyng et al. (1981) reported images in the 
energy range 16-30 keV which change from a double footpomt structure early in 
the event to a more amorphous structure later. The Japanese free flyer, YOHKOH, 
is presently making observations; however, the data have not been reduced yet. 

Figure 1.12 shows HXIS hard x-ray images of a flare on May 21, 1980 (An- 
tonucci et al., 1982), which are somewhat ambiguous, due to the limited temporal 
and spatial resolution of the instrument among other reasons. However, these and 
other HXIS hard x-ray images indicating the existence of footpoints as sources do 
not sufficiently discriminate between the Model T and the Model N-T. In fact, due 
to low count rates and the correspondingly large integration times necessary to 
achieve a statistically significant image, and to telemetry limitations (MacKinnon, 
Brown, and Hayward, 1985), the HXIS data has in general done little to alleviate 
the uncertainty between the two models (Brown and Emslie, 1988; Emslie, 1983). 
More statistically significant hard x-ray imaging data, at significantly higher photon 
energies and with higher resolution, are vitally needed to provide better descriptions 

of these processes. 

The spatial resolution required to understand the physical processes at work 
in the flare may be calculated as follows. From Emslie and Machado (1987), the 
energy of the incident electron may be related to the column depth to which it 
penetrates a thick target by 


E 2 = El - 2KN C 1 - 1 ) 

where E is the electron energy at particular column depth, N = fndz, K = 
2.8 x 10 -18 cm 2 keV -1 , and n is the electron density. We can estimate the depth 
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Figure 1.12: HXIS Hard X-ray Images 

This illustrates images of a flare in hard x-rays taken on May 21, 1980 (Antonucci 
et al., 1982). Although the resolution is somewhat limited, footpoints seem to be 
indicated. 
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an electron will penetrate into the chromosphere from 


El — E 2 

r sj O 

z = 


2Kn 


( 1 . 2 ) 


Taking E a = 20 keV and E = 10 keV (the threshold for production of a 10 
keV photon), we find with n = 10 n cm -3 , 


^ _ x = 7 arc second. (1-3) 

10 11 

Therefore, in order to follow the degradation of the electron energy, a spatial so- 
lution substantially better than 7 arc seconds is required. Hence, the minimum 
requirements suggested by this analysis for imaging telescopes are spatial resolu- 
tions on the order of 1-4 arc second. 

Sub arc second resolution may be desirable in the next generation of instru- 
ments to further characterize these processes. Similarly, from considerations of 
collisional stopping times, and loop transit times for a typical loop 1 arc second in 
cross sectional diameter and 20 arc seconds in loop diameter, temporal resolutions 
on the order of 1-5 seconds are desirable (LaRosa, 1990). 


C. X-Ray Imaging Background 

A Fourier-transform telescope typically consists of a set of subcollimators each 
of which modulate incoming photons onto an associated detector array providing a 
measurement of a Fourier component. Image making in hard x-rays can be shown 
(Chapter IV) to correspond to that in the radio frequency regime in that the data 
are basically now in the same form as those from a radio telescope. Also, it can be 
brought into the Astronomical Image Processing System (AIPS), hitherto used for 
radio images and processed with the algorithms therein. Refer to Chapter IV for a 
more detailed discussion. 
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AIPS, developed and maintained by the National Radio Astronomy Observa- 
tory (NRAO), is widely accepted today in the radio astronomy community (Mertz, 
1989). It is a highly sophisticated VMS overlay software system which requires 
over 180,000 blocks of VAX disk space to rim. At its core is a multitude of software 
routines which allows the user to bring Flexible Image Transport System (FITS) 
data into the AIPS operating system; (Greisen, 1986; Greisen and Harten, 1981) 
fast Founer transform the data to form a dirty map; clean the dirty map to form 
an image using several different, optional algorithms; and then display the image in 
a number of different ways. Peripheral software offers the user a number of options 
with which to enhance an image in a number of different ways. Chapter IV will 
discuss this in greater detail. 

One major advantage in using AIPS is that it allows a “dirty” image to be 
cleaned several different ways. Two of the most widely used are the CLEAN algo- 
rithm and the Maximum Entropy Method (MEM). Figure 1.13 shows an example 
of a dirty image of an astrophysical source emitting in the radiofrequency regime. 
Although many artifacts are present, the trained observer will expect a real source 
or sources at the center of the symetrical pattern. The shape and intensity of the 
source is basically unknown. Figure 1.14 shows the true image after the CLEAN 
algorithm was applied. Here, the source, a radio emitting galaxy, is revealed much 
more clearly. 

One of the strengths of a numerical model is that a dirty image can be cleaned 
conveniently using both techniques and then the images compared. In addition, ar- 
tifacts which might be introduced by approximations in the techniques will become 
apparent and will serve as useful supporting information for handling real data. 

D. Research Approach 

We used the wealth of literature relating to both Model T and Model N-T 
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Figure 1.13: A Dirty Image Of A Typical Astrophysical Source 

The trained observer will expect to find a source at the center of the symmetrical 
region. However, this dirty image alone is not sufficient to show accurately the 
shape and intensity of the source. 
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Figure 1.14: The Clean Image 


The CLEAN algorithm was applied to generate this cleaned image. When 
paring the quality of this image with its dirty image, the power of this i 
reconstruction algorithm becomes truly evident 



to describe quantitatively the hard x-ray emission as a function of time, energy, 
and loop position (Chapters II and III). The emission profiles for each model are 
then examined for distinctive signatures and used to construct synthetic flares to 
be viewed by the model telescopes. 

Using random photon counting techniques, SMC and RMC models were con- 
structed to view Model T and Model N-T flare profiles and AIPS was modified and 
calibrated to reconstruct images from the telescope x-ray optics output. Compari- 
son of the images at similar times and energies should provide a basis for judging 
the ability of a given telescope to resolve the distinctive features of a synthetic 

event. 

In addition, various instrumental questions were addressed such as the effects 
of low signal levels, random noise on the detector, twist, etc. The end-to-end set of 
software including those codes which generated synthetic flares, those codes which 
represented the operation of the telescope, those codes which developed the appro- 
priate detector responses, and those codes in AIPS were designated the Marshall 
Space Flight Center Hard X-ray Imaging Simulations (MAHXIS). See Appendices 
A - D for examples. 
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II. MODEL N-T HARD X-RAY EMISSION 


A. INTRODUCTION 

As mentioned in the introduction, there axe two prevalent, competing theories 
(Model T and Model N-T) for explaining hard x-ray emission from solar flares (e.g., 
Campbell et ah, 1990; Brown, 1974; Emslie and Vlahos, 1980; Tandberg-Hanssen 
and Emslie, 1988). The purpose of this chapter is to describe quantitatively the 
Model N-T hard x-ray emission as a function of time, energy, and loop position. 
This will be compared with the Model T profiles to be discussed in Chapter III 
to determine if any differences exist between the two. These profiles will then be 
viewed numerically (Chapter VII) through a Fourier telescope model to determine 
if in fact the differences will be evident in the hard x-ray images. 

Figure 2.1 illustrates the Model N-T concept. A beam of electrons is injected 
at the apex of the loop. These electrons travel down the legs of the loop and impact 
the denser chromosphere in the footpoints. Since the bulk of the interactions in 
this region are Coulomb collisions, the plasma is heated rapidly which increases the 
pressure locally. This increase in pressure forces material to migrate up the legs of 
the loop until it reaches the apex after which point it simply piles up in the loop, 
increasing the density. This migration leads to substantial increases in temperature 
and in density in the loop. As the density increases, more and more interactions 
between beam electrons and the ambient plasma cause additional heating and also 
result in nonthermal bremsstrahlung emission. The hot plasma in the loop eventu- 
ally reaches a temperature in which thermal emission becomes significant at hard 
x-ray energies. 

The Model N-T exhibits both thermal and nonthermal characteristics. The 
distinction between thermal and nonthermal bremsstrahlung used in this context 
is as follows: if the emission is produced by interactions within a population of 
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Figure 2.1: Model N-T Concept 

A beam of electrons is injected at the apex of the loop (graphics from Tandberg- 
Hanssen and Emslie, 1988). These electrons travel down the legs of the loop and 
impact the denser chromosphere in the footpoints. Since the bulk of the interactions 
in this region are Coulomb collisions, the plasma is heated rapidly which increases 
the pressure locally. This increase in pressure results in the migration of material 
up and down the legs of the loop. Both nonthermal and thermal emission occur 
throughout the loop at lower energies. 
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particles roughly distributed according to a Maxwellian then we say the process is 
thermal in nature. On the other hand, if the population forms a separate high- 
energy component at energies much greater than the bulk of the population then 
we define it as nonthermal. 

Brown and Emslie (1987) calculated the spatial structure and temporal evolu- 
tion of the hard x-ray emission of both nonthermal and thermal components during 
the impulsive phase of an electron- heated (Model N-T) solar flare based on hydro- 
dynamic simulations by Nagai and Emslie (1984) and Emslie and Nagai (1985). 

Recently, Mariska, Emslie, and Li (1989) carried out a series of numerical 
simulations which more completely describe the response of the solar atmosphere 
to heating by a nonthermal electron beam. In this work, the electron beam is 

injected at the apex of the loop and the beam electrons have an energy spectrum 
of the form 


F 0 (E 0 ,t) = Fo(E 0 )g(t) 

(2.1) 

4(6-2)F ig (t)[E 0 ] 2 
(S + 2 )£? UJ 

(2.2) 

E 0 55 E c 



or 


if 


F 0 (E 0 ,f) 


4QS-2 )F,g(t) pel-* 
(S + 2)EI [E c 


(2.3) 


E 0 > E c . 
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Here Ft is the injected flux (erg cm ~ 2 s l ),E c is a low-energy inflection point , 6 is 
the spectral index, E 0 is the electron energy in keV, and g(t) is a triangular time 

profile given by 



(2.4) 


for 


0 < t < 30 s 


or 


»« = 2 -To 


(2.5) 


for 


30 < t < 60 s. 

The low-energy form of the electron energy spectrum (equation 2.2) gives a 
much smoother variation of heat input with depth than that corresponding to an 
electron spectrum with a sharp low-energy cutoff (see e.g., Nagai and Emslie, 1984) 
and is primarily used to permit a more efficient and rapid numerical solution of the 
hydrodynamic equations. Its use does not change, however, the power law spectral 
shape of hard x-ray bursts at photon energies (e > E c ) and thus is still consistent 
with observations (Figure 1.10). 

Figure 2.2 shows the results for the density in the model from 10 to 30 seconds. 
At 10 seconds into the flare, bulk material motion can be observed moving up from 
the footpoints. At 20 seconds the material has almost reached the apex of the loop, 
and at 30 seconds the material is piling up in the legs of the loop. The change in 
density from 20 to 30 seconds produces some interesting spatial structure for the 
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emission at 10 keV (Figure 2.6) since both the thermal and nonthermal emission 
are directly affected by changes in density at this energy which will be discussed 
later. 

Now, following Emslie and Machado (1987), the integrated nonthermal hard 
x-ray yield from column depth 0 to N (N= / ndz) in photons cm -2 s -1 keV -1 may 
be written as 


where 


f N 

J(e,N)= / I(e,N')dN' 

Jo 


6 Jn '= o Jez 


Fp(E o ) 

S {El - 2KN ) 


ycIEq 


( 2 . 6 ) 


(2.7) 


E; = (e 2 +2KN)i 


( 2 . 8 ) 


emd 


K = 27re 4 A (2.9) 

is the coefficient of the variation of energy with column depth due to Coulomb 
collisions and is defined in Tandberg-Hanssen and Emslie (1988). A is the Coulomb 
logarithm. From Brown (1971) and Li and Emslie (1987), the cumulative intensity 
J(e,N ) may be written as 

(<+2 -';(*; v) (2.io) 

in umts of photons cm -2 s — 1 keV~ 1 at the Sun, where 

A = 4 x 1.6 x lO -9 -^: = 1140 erg* 1 . (2.11) 
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Figure 2.2: Density Versus Loop Position 

This shows the density results from the Mariska, Emslie, end Li (1989) hydrody- 
namic simulation as a function of loop position for a family of times. The beam 
heats the material in the footpoints and the bulk motion from the footpoints to- 
ward the top of the loop can be seen. At 30 seconds, the front has reached the top 
and material starts piling up. This causes approximately an order of magnitude 
increase in the density from 20 to 30 seconds at the apex. 
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The dimensionless quantities 


2 KN 


( 2 . 12 ) 


and j are functions of both the dimensionless column density x and y = -=r. These 

Cj c 

axe shown graphically for various forms of F 0 (E 0 ) by Li and Emslie (1990), and 
Emslie and Machado (1986). 

As a result of Coulomb collisions of the beam electrons with the ambient 
electrons in the target, the electron temperature in the flare loop can rise to 
100 x 10 K (Mariska, Emslie, and Li, 1989), adequate to produce signifi- 
cant hard x-ray emission, especially at lower energies (e < 10 keV). This additional 
secondary” thermal bremsstrahlung yield from a segment of the loop is (Allen, 
1973; Culhane and Acton, 1970) 


A j w = 12 * 10 11 (^f )"*(*?) < 213 > 

in photons cm 2 s 1 keV 1 at the Sun, n, T , and l are the density, temperature, 
and the length of the segment. 

A given AJ (whether due to nonthermal or thermal emission) produces a cer- 
tain number of counts on a detector in low-Earth orbit can be obtained from 


C = Ct A J (counts s *), (2-14) 

where 


CT = a( 4^ )At ’ < 2 ' 15 ) 

and A is the flare area, R = 1 AU, o is the effective detector area, and Ae is the 
effective energy width of the observation bin. 
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This equation assumes an ideal detector which is equally responsive to photons 
of any energy within the selected bin. In reality, above 100 keV detector quantum 
efficiency tends to decrease dramatically as one goes higher in energy; thus, it is 
desirable to take energy bins as small as possible to compensate for this effect. Also, 
the detector response characteristic may change over time. Real-world detector 
contracts will be discussed in greater detail in Chapter IV. 

In the following section the total hard x-ray bremsstrahlung emission as a 
function of energy, time, and position along the loop is presented. The beam 
parameters have the following values: F t = 10 11 erg cm -2 s \ E c = 15 keV, and 
6 = 6 . 

B. EMISSION AS A FUNCTION OF TIME, POSITION, AND ENERGY 

The following series of plots illustrate the predicted emission from model N- 
T at the Earth. This, as discussed previously, will be composed of both thermal 
and nonthermal components. However, since our telescopes cannot, of course, dis- 
tinguish how an x-ray photon was created we simply look at the total emission 
signature. The following plots will show the intensity at the Earth in photons 

cm -2 sec -1 keV -1 for a given loop pixel. 

The first two plots (Figures 2.3 and 2.4) show the predicted spatial distribution 

of the emission due to thermal effects only for a given family of times (10, 20, 30, 
40, 60 seconds), for photon energies at 10 and 40 keV, respectively. As one would 
expect, the decrease in thermal emission from 10 to 40 keV is extremely large. 
At 40 keV and higher there is no significant emission being produced by thermal 
bremsstrahlung. This is to be expected from the relatively low temperatures (kT 
< 10 keV) found in the loop at this time (equation 2.13). 

Figure 2.5 shows the predicted emission profile of the flare at 10 keV due to 
nonthermal emission alone for the same family of times mentioned above. While 
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Figure 2.3: Thermal Emission Profile at 10 keV 


This shows the Model N-T thermal emission component only. The thermal emission 
starts increasing between 20 and 30 seconds, corresponding to the arrival of bulk 
material from the foot points, and remains high throughout the rest of the event. 
At 60 seconds, the beam has shut off and only the emission from the filled loop 


remains. 
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Figure 2.4: Thermal Emission Profile at 40 keV 

The maximum temperatures in the loop (=2x 10 7 K) are not high enough to 
produce significant thermal emission at 40 keV. Hence nearly all the significant 
emission from Model N-T at energies from 40-100 keV is nonthermal. 
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most of the emission occurs in the footpoints, significant emission is also produced 
higher in the loop. Corresponding to the large density increase at the top of the loop 
from 20 to 30 seconds is a substantial increase in the total emission from the apex 
of the loop. The beam heats this material producing some nonthermal emission 
while an order of magnitude additional thermal emission is produced by the heated 
plasma. This is to be expected since the thermal bremsstrahlung is i nh erently a 
more efficient process (Chapter I). At 60 seconds (beam shutoff), only a filled loop 
remains. 

Figures 2.6 through 2.9 show the combined emission at the aforementioned 
times for 10, 40, 70, and 100 keV photons, respectively. The footpoint region 
is heated rapidly by the beam and the corresponding increase in local pressure 
produces both upward and downward mass motions. This appears in Figure 2.2 as 
a shift in the high density region downward and is similarly reflected in the emission 
by a corresponding shift in the peak intensity from one pixel to the next lower one. 

As stated before, the emission, seen in Figures 2. 7-2.9, at 40-100 keV is non- 
thermal; and the profile signature is one of relatively bright footpoints with less 
bright pixels in the upper part of the flare loop. 

C. SUMMARY 

The emission as a function of time, energy, and loop position of the Model 
N-T has been determined quantitatively. The predominant hard x-ray signature 
predicted by this model is one of twin footpoints early in time at low energies 
and at all times at higher energies. The telescope will image both thermal and 
nonthermal components together; thus, the total profile must be used as input 
for the Fourier telescopes, the SMC and the RMC, discussed in Chapters V and 
VI. How well these telescopes can image the predicted profiles will be discussed in 
Chapter VII. 
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Figure 2.5: Nonthermal Emission Profile at 10 keV 

Clearly, emission from the footpoints is dominant in the loop from 0-20 seconds and 
remains at significant levels until beam cutoff at 60 seconds. However, significant 
emission is also produced higher in the loop especially after those times m which 
the density has greatly increased. 
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Figure 2.7: Total Emission Profile at 40 keV 

This shows the Model N-T total emission profile. These combined thermal and 
nonthermal profiles will serve as input to the telescope models (Chapter VII). The 
movement in time of the footpoint peaks corresponds to mass motions both up 
and down the loop due to increased local pressure in the region. As expected, the 
total intensity is significantly reduced at 40 keV relative to that at 10 keV and is 
predominantly nonthermal. 
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Figure 2.9: Total Emission Profile at 100 keV 

This shows the Model N-T total emission profile. As expected, the total intensity 
is significantly reduced at 100 keV relative to that at 70 keV (Figure 2.8) and is 
predominantly nonthermal. 
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Clearly, the Model N-T evolves rapidly over a period of a few seconds. The 
bright footpoints spatially occupy 1-3 arc second each and should be observable 
by a telescope with 2 4 arc second resolution such as ours. While the footpoints 
persist in this model for several seconds, having 1 second integration times would 
still be desirable as the 10 keV case evolves significantly from a bright footpoint 
structure to a predominant central peak and then to a filled loop. 
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III. MODEL T HARD X-RAY EMISSION 


A. INTRODUCTION 

Chapter II discussed the quantitative development of the Model N-T hard 
x-ray emission as a function of time, loop position, and energy. This chapter will 
present the Model T concept and determine analytically its associated emission 
characteristics. Once the emission profiles for this model are obtained, they can 
be compared with those in Chapter II to determine if those for each of the two 
models are sufficiently distinct to afford an observational test. Chapter VII will 
discuss Fourier telescope image differences between the two models. 

The Model T concept is one of an extremely hot mass of plasma created in 
a coronal loop such that the thermal electrons have sufficient energy to produce 
bremsstrahlung at hard x-ray energies through large-angle scattering off ambient 
protons (Chubb et al., 1966). Based in part on fusion research, Brown et al. 
(1979) describe an improvement to this basic idea. Initially, the electrons m a 
small region of the loop near its apex are heated rapidly. The heating mecha- 
nism and rate are assumed to be such that only the electrons would be heated 
appreciably leaving the ions at essentially their preflare (= 10 6 K) temperature. 

A large electrical current is generated by the hot electrons streaming out- 
ward from the relatively small, initial source volume down the legs of the flare 
loop containing cooler plasma. The electric field established by a combination 
of inductive and charge separation effects accelerates ambient electrons into a 
reverse current which quickly exceeds the ion-sound velocity in the plasma, thus 

generating ion-acoustic waves (Figure 3.1). 

A brief discussion of the physics of ion-acoustic waves follows (Krall and 
Trivelpiece, 1973; Chen, 1974). In the absence of ordinary collisions within a 
plasma, ions can still transmit momentum to each other through the intermediary 
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Figure 3.1: Model T Concept 

The electrons in a relatively small region at the apex of the loop axe heated 
to 2 x 10 8 K. The hot electrons begin streaming away down the legs of the 
loop inducing a return current exceeding the ion-acoustic speed which leads to 
ion-acoustic turbulence. The waves form a marginally stable, moving barrier 
which contains most of hot electrons. The most energetic ones escape and serve 
to maintain the marginal stability condition. Thermal emission is produced in 
the kernel while nonthermal emission is produced outside the kernel and in the 
footpoints. 
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of an electric field carried by mobile electrons. These form regions of compression 
and rarefaction, just as in an ordinary sound wave. Compressed regions tend to 
expand into the rarefactions. This happens because (i) thermal effects tend to 
spread out the ions and (ii) the ion bunches are positively charged and tend to 
disperse because of the resulting electric field. Once in motion, the ions in the 
bunches tend to overshoot equilibrium because of their relatively large inertia. 
This naturally leads to an osciUation in which the compressions and rarefactions 
axe regenerated to form an ion-acoustic wave. 

Ion-acoustic waves occur due to the interaction of both the electron and ion 
components within the plasma. These waves are strongly affected by Landau 
damping and hence only created in any strength when T e > T t . If one considers 
the classic surfboard analogy (Chen, 1974) in which surfboards (electrons) slower 
than the (ion) wave gain energy from the wave while surfboards faster than the 
wave lose energy to the wave, then the importance of the electron temperature 
being higher than that of the ions becomes evident as this spreads out the electron 
distribution and damps the wave more slowly. 

Figure 3.2 illustrates the two distributions at the beginning of the event near 
the top of the loop in the Model T. The hot electron distribution is skewed repre- 
senting a return current moving through the ambient, colder ions. The threshold 
for instability is reached when the bulk velocity of the backward- streaming elec- 
trons exceeds the ion-acoustic speed. 

Rapid isotropization of ion-sound waves occurs through scattering on the 
ambient particles, and the resulting ion-sound turbulence in turn scatters the hot 
electrons. This turbulence in effect forms a moving barrier in which the electrons 
have a reduced mean free path. Thus, most of the hot electrons in the kernel are 
reflected while only the ones possessing the highest energies (i.e., greatest mean 
free paths) pass through. This passage of the highest energy electrons is vital to 
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Figure 3.2: Electron and Ion Distribution Functions 

This shows the cold ion and hot electron distributions near the top of the loop 
at the beginning of the event in the Model T concept. The electron distribution 
at low velocities is skewed reflecting a return current. When the bulk velocity of 
this return current exceeds the ion-acoustic speed, an ion-acoustic instability is 
excited leading to ion-acoustic turbulence. 
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maintaining the marginal stability of the barrier throughout its movement from 
the apex to the footpoints corresponding to the flare’s impulsive phase. 

Marginal stability in this context may be illustrated by considering two 
perturbations. In the first case, the level of acoustic wave turbulence is increased 
slightly. This causes the hot plasma to be contained more efficiently reducing 
the escape current and the return current and hence the level of acoustic wave 
turbulence. In the second case, the level of turbulence is decreased slightly. This 
reduces the effectiveness of containment and allows the current to increase. The 
return current then increases restoring the level of turbulence to its original value. 

A simple analog of marginal stability is presented in Figure 3.3. In this ana- 
log, a system comprised of hose pipes, movable barrier, and bucket are arranged 
as shown in the figure. The hose pipes have their nozzles adjusted so that the 
water flow rate through each one varies. The system is placed in equilibrium 
by placing the barrier so that the fastest streams of water clear the top of the 
barrier to be collected in the bucket. A drain hole is punched in the bottom of 
the bucket. Some water clears the barrier and is collected by the bucket adding 
to its weight. Thus, the barrier is maintained at its equilibrium height balanced 
between the pull of gravity and the spring force. 

If the barrier is pushed down slightly, more water goes to the bucket, which 
gets heavier and pulls the barrier back up again. If the barrier is pushed upward 
slightly, less water gets over the top and the bucket gets lighter thus allowing the 
spring to pull the barrier back down again. Finally, we note that if the fastest 
hose is turned off, then the bucket gets lighter allowing the spring to pull the 
barrier down slightly to allow additional water over the top thus maintaining 

equilibrium. 

Returning to the solar case, two barriers would form, one in each leg of 
the loop and would move apart down the legs of the loop at about the lon- 
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Figure 3.3: Marginal Stability Analog 

This shows streams of water projected at a barrier supported by a bucket and 
pulley system. Some of the water gets over the top of the barrier into the bucket 
but flows out through a hole in the bottom. The system can be designed to be 
initially in equilibrium. If a high velocity stream is turned off, the barrier will 
drop down a bit to maintain equilibrium. 
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sound speed. The region in between would contain a hot kernel of expanding 
plasma. Brown et al. (1979) stated that these barriers would be relatively thin 
(on the order of 0.13 km). Other authors (Smith and Brown, 1980) have revised 
this estimate upward somewhat; nevertheless, all agree that the thickness of 
the barriers is relatively small compared to the overall length of the loop and 
therefore not resolvable. 

Electron and ion densities would remain roughly constant over the loop to 
within 10-20 % (Smith and Lilliquest, 1979). This may visualized simplistically 
as a weather front in our own atmosphere. There is relatively little change m the 
total density through the front, yet, as it moves, the temperature at a particular 
location may change significantly. The predominant loss of energy would be 
in the form of an anomalous heat flux carried by the escaping tail electrons, 
Q -= n e m e vlc a . 


B. EMISSION CHARACTERISTICS 

The aim of this chapter is to describe the Model T hard x-ray emission. Once 
the emission in time and energy is known, one can then turn to the telescope 
simulations to determine at which times and at which energies this hard x-ray 
model would provide a discemable and characteristic signature to the observer. 

While it is not necessary to rework here the details of the internal processes 
which form and maintain the barrier, it is necessary to note that it is permeable 
first to the highest energy available electrons. This happens because only the 
highest energy particles have sufficient energy (i.e., long enough mean-free-paths) 
to diffuse all the way through the barrier region (Brown et al., 1979). Another 
important implication is that a collisionless flux stream which is not energetic 
enough to escape early in time will have another chance later simply because it 
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would then represent the highest available energy. In the meantime, however, it 
is reflected by the barrier. 

Within the kernel, hard x-rays are produced by thermal bremsstrahlung. In 
the footpoints and to some extent in the legs of the loop, a thick target, non- 
thermal description is more appropriate since we have a distribution of electrons 
with energies E >> kT impacting the relatively cool, lower corona and chromo- 
sphere. While there are many initial conditions that could be chosen to quantify 
the treatment, the ones used here form a typical case. 

The Model T basic description is restated as follows. The initial electron 
distribution is a Maxwellian at 2 x 10 8 K. The initial kernel length is 10 8 cm 
bounded by relatively thin regions of confining ion-acoustic turbulence. The 
centroid of the initial kernel is taken to occur at the apex of the loop. An 
unknown mechanism heats the electrons initially and is assumed to have little 
effect on the more massive ions (i.e., T e T t ). Electron and ion number densities 
are taken as constant at 10 11 cm -3 throughout the impulsive phase (Smith and 
Lilliequest, 1979). The loop geometry is a semi-torus with a = 1 arc second 
cross section and an overall loop diameter of = 20 arc second. The speed of the 
turbulent acoustic wavefronts is taken to be the local ion-acoustic speed (Brown 
et al., 1979). 

Given these assumptions, we are now in a position to determine T e , the 
temperature of electrons inside the kernel, and L, the length of the kernel as 
functions of time. Assuming energy losses from the kernel to be relatively small, 
conservation of energy may be used to write 


LT e 2* L a T n 


The source length L increases at the local sound speed, c s . 


(3.1) 
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(3.2) 


e .1 

C s = ( ) 2 

mi 


on each side; thus 



which by (3.1) may be written as 


dL_ 2 j kLpf^ 

dt V m iL 

Integrating, we obtain 


(3.3) 


(3.4) 


= L 0 

3 

1 + \ 

/“■• ; 


Lq \ 

1 mi . 


(3.5) 


Using equation (3.1) again, the electron temperature as a function of time may 
now be written as 


T e (t) = T 0 


1 



(3.6) 


Numerically, for L 0 = 10 8 cm and T 0 = 2 x 10 8 K, the time for the barrier 
to reach the footpoints (i.e., the length of the impulsive phase) is 28.5 seconds 
which is consistent with typical observations for the duration of the impulsive 
phase (Tandberg-Hanssen and Emslie, 1988). The final temperature is 8.6 X 10 6 
K is sufficiently large compared to the preflare temperature (on the order of 10 6 
K) that equation (3.1) is a valid approximation (Figure 3.4). 

The aforementioned discussion is concerned with the bulk behavior of the 
distribution; however, the high-energy tail of the electron distribution must also 
be taken into account as it will be responsible for any emission in the footpoints 
of the loop. From Spitzer (1962), the collision time, t c , may be written as 
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Figure 3.4: Model T Kernel Temperature As A Function Of Time 

This shows the evolution of the bulk temperature of the expanding kernel as a 
function of time. It cools very rapidly early in the event and then slows expo- 
nentially for the later portion of the event. One would expect from this behavior 
that the high-energy emission will come early in time. 
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(3.7) 


_ ,266Tl 

n/nA ’ 

where In A is the Coulomb logarithm, which for solar conditions may be taken 
as a constant (In A = 20; Spitzer, 1962). Since the self collision time, t c , for the 
electrons is short (« 10" 1 s) relative to the characteristic expansion and cooling 
time of the kernel (« 30 s), one expects the kernel for the bulk of the electron 
distribution to be Maxwellian (Figure 3.5). However, early in the impulsive phase 
one would expect the tail of the distribution to be non- Maxwellian (Kahler, 
1971a, b). Since the initial Maxwellian tail is essentially collisionless over the 
length of the kernel at that time, one would expect it to persist for some finite 
time while the bulk of the distribution follows a relaxing Maxwellian throughout 

the event. 

In this regard, MacDonald et al. (1957) discuss the behavior of a system of 
particles relaxing through Coulomb forces from an extreme initial distribution 
to a final Maxwellian. By treating the problem with dimensionless parameters 
(e.g., dimensionless velocity £)> the results are applicable to a distribution at any 
temperature. As a specific example, they chose an initial distribution in the form 
of a displaced Gaussian and provided a few curves at different times, including 
in particular the one for the distribution after a relaxation time of 5.9f c . At this 
time, the relaxing distribution approximated a Maxwellian except in the tail. 
The temperature used to calculate the collision time was that associated with 
the fully relaxed distribution. 

This specific example may be generalized to be applicable to our situation 
by defining the dimensionless differences 


D 0 id = 


hi — h 

hi — h\f 


(3.8) 


and 
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Figure 3.5: Maxwellian Cooling Over Time 

The kernel cools slowly enough that the bulk of the distribution will follow the 
evolution of a cooling Maxwellian distribution. However, the tail of the distribu- 
tion will not necessarily be Maxwellian. 
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h — h,M 
hi - hM 


(3.9) 


representing the difference between the initial, final, and Maxwellian distributions 
as a function of velocity (Dm + £>„«„ = 1). Here h is the distribution after a 
relaxation time of 5.9 t c , ku is the fully relaxed Maxwellian, and h, is the initial 
distribution (Figure 3.6). 

Except in the region where hi « h M at £ = 0.25 (Figure 3.6), D old and 
D new are expected to apply to a more general class of problem such as ours. 
If we plot D old and D new (Figure 3.7) against the dimensionless velocity, we 
find that D new = 0 over the bulk of the distribution and then begins to rapidly 
grow as one approaches the tail of the distribution. Unfortunately, not enough 
data are provided to determine a similar departure from the initial distribution. 
Now, a velocity may be selected ( D ncw = 0.01) as the cutoff at which the relaxing 
distribution is defined to be no longer Maxwellian. This point on the MacDonald 
et al. (1957) curves corresponds to 


6 = 0.31 


(3.10) 


at time 


t = 5.9 t c ( 3 - n )’ 

where 

f = 0.335 (3-12) 

is the dimensionless velocity used by MacDonald et al. (1957), T e is the tempera- 
ture associated with the new distribution, k is the Boltzmann constant, and m is 
the electron mass. Equation (3.12) can be inverted to provide the velocity in cm 
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Dimensionless velocity £-0.335 v'/mpTfcT 


Figure 3.6: Relaxation To Maxwellian Distribution 

This shows the initial, dimensionless Gaussian distribution and four later times 
as it relaxes toward a final Maxwellian. Note that for the dimensionless time 
484.17 (i.e., 5.9 t c ), the curve closely approximates that of a Maxwellian for 
lower values of f (MacDonald et al., 1957). 
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s- 1 which defines the edge of the Maxwellian bulk distribution in phase space, 
Wl . We also have sufficient information at this point to derive an expression for 
the temperature of the bulk distribution corresponding to 5.9 t c and hence u e at 

that time. 

From equation (3.6), 



(3.13) 


Setting t = 5.9 t c , the time used by McDonald et al. (1957), we find, using 
equation (3.7) and the substitution 


X 



(3.14) 


that equation (3.13) may be written in the form of a quadratic equation 


ax 2 + bx + c = 0 


(3.15) 


where 


a = 5.9 x 


0.266 
nln A 


(3.16) 


_ I / m < r 

3 y LokTo ° 


(3.17) 


and 


-u 


mi 

L o kT 0 


(L 0 T 0 ) 2 . 


(3.18) 


Hence, equation (3.15) can be solved to yield T e at the time corresponding 
to the MacDonald et al. (1957) results on which Figure 3.6 is based. This also 
will allow a determination of the rms velocity, 
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Figure 3.7j Dimensionless Differences 

This shows a comparison of the dimensionless differences, D new and D oU ob- 
tained from the curves provided by McDonald et al. (1957). The point at which 
Dnew exceeds 0.01 was chosen as the cutoff for the bulk Maxwellian. This cor- 
responds to a dimensionless velocity of £ = 0.31 at 5.9 t c . Clearly, above this 
cutoff, Dnew begins to diverge from a Maxwellian rapidly. 
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(3.19) 



The dependence of collision frequency on velocity may be used to describe fully 
the electron distribution as a function of time throughout the event. We define 
Vcrit as the velocity at which the number of collisions occurring during a partic- 
ular elapsed time corresponds to 1 (N(v, t ) = 1). Spitzer (1962) decribes one 
collision as the summed effect of many grazing encounters by a particular parti- 
cle until a 90° change in direction has been accomplished. Clearly, at velocities 
less than Wl , defined by equations (3.10) and (3.12), enough of these collisions 
have occurred for the bulk distribution to be Maxwellian. We therefore expect 
this velocity to satisfy the relationship vi < v crit . Similarly, at a higher velocity, 
V2 > v crit , the distribution must resemble the initial distribution as beyond that 
point collisions are not occurring quickly enough to change the character of the 

distribution in this region of phase space. 

Now, the number of collisions at a given velocity, N, is given by 


jV(v, t) = vt ( 3 - 20 ) 

where t is the time and u( v) is the collision frequency associated with a given 
velocity. The collision frequency varies with velocity as 

„ = ^ (3.21) 

t ; 3 

where i/ e = is the collision frequency associated with the rms velocity, iv 
Substituting for v in equation (3.20) we find 



(3.22) 
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Substituting the value of v x = [from equations (3.10) and (3.12)] 

for v, and the value of t = 5.9t c from equation (3.11) provides the number of 
collisions required to define the boundary ( D new = 0.01) in phase space of the 
Maxwellian bulk distribution. This is 


Ni = 5.9 



3 


s 1.4, 


(3.23) 


and indicates that only a relatively small number of Spitzer collisions axe needed 
to thermalize an initially non-Maxwellian bulk distribution. Since uj £ (equa- 
tion 3.22), and Vi = ( 0 . 335 * ) Ve t = 5.9t c , we obtain an expression for as a 
function of time; 


Substituting we find that 





(3.24) 


”1 = (7 x 10 9 )*3. (3.25) 

This finding agrees with MacDonald et al. (1957) who found that the rate of 
expansion of the edge of the Maxwellian bulk distribution in phase space is 
independent of temperature. 

Now, for N = 1, v = v cr it , and we have 


t_ 

U 



At t = 5.9 t c , we therefore find that 


(3.26) 


Vcrit = (5.9)3 U e ^ 1.8 V e 
and from equation (3.24) at the same time, 


(3.27) 
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Vi = 0.9(5. 9)3v e = 1.6 v e 


(3.28) 


Now, in the absence of information from MacDonald et al. (1957) for the 
high energy part of the distribution, the number of collisions, N 2 , beyond which 
the distribution remains that of the initial condition is not well known. However, 
we observe that tq < v crxt < v 2 and that by assuming symmetry in the difference 
curves shown in Figure 3.7, we derive 


2 Vcrit =V X +V 2 ( 3 ‘ 29 ) 

and 

N 2 = [(2AU,)* - - °' 7 ' (3 ' 30) 

Also, the corresponding velocity at 5.9 t c will be 

v 2 = l.l(5.9)3 Ve s 2.0 v e . ( 3 - 31 ) 

Thus, for velocities, t> < «*, the distribution is a relaxed Maxwellian at the 
new T and for v > v 2 the distribution is represented by the initial condition. The 
middle region is relatively small in extent and may be interpolated by a cubic 
spline which matches endpoints and first and second derivatives. 

Previous authors have assumed all escaping electrons to disappear from the 
kernel instantaneously. Clearly, however, since L = 10 8 cm, these collisionless 
electrons will not all arrive at the front instantaneously but rather will replenish 
the front electrons continuously during some finite time interval. We will now 
treat this replenishment more correctly in order to better understand its role in 
the event. We start with the anomalous heat flux equation from Brown et al. 
(1979) to calculate the lower cutoff velocity of the escaping electrons, v lco , at a 
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given time. Knowing the velocity as a function of time is the key to describing 
the escaping flux traveling to the footpoints. We write the heat flux as 


fVmax 

Q= F(v,t)v 3 dv 

Jvico 


where F(v,t) is the electron distribution. 
Now, for marginal stability, 


(3.32) 


Hence, 


Qms = 


3 

-nm e v 


2 

e 


Cs 


(3.33) 


3 

-nm e v e 2 c s = / v z F(v,t)dv 

'Vico 

where both v max and vi co are functions of time. 

For the initial condition, we have v max (t) = v maxo 
first vico we find that 


(3.34) 


oo. Solving for the 


where 




c 3 m e 


3 (Wjf 

2 m.i 2 


f(v, 0 ) = 



Evaluating the integral I, we obtain, 


(3.35) 


(3.36) 


1 2 G^t) eXp { [(If) + (^7)" 

a 

Setting this result equal to | ^-^7 2 leads to the equation 

m 2 


2 

Ico 


(3.37) 
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(3.38) 


c 2 + v, co 2 = 3ci exp 



where Ci and c 2 3-re constants; 


and 



Solving, we find that 


2kT 

c 2 = • 

m e 


(3.39) 


(3.40) 


v,„ = 2.6t)«, ( 3 - 41 ) 

in agreement with Brown et al. (1979). 

While this works for the Maxwellian initial condition, we would not expect 
it to be valid for our distribution at later times since v max will become finite as 
the highest energy electrons escape. For later times consider the following: the 
collisionless flux may be modeled as a series of differential flux streams x, with 
velocities between v, and v t + dv. These streams are continuously striking the 
front and either (a) passing through to maintain the marginally stable turbulence 
(if u > v lco ) or (b) being reflected if v < v, co . As soon as the highest velocity 
stream departs, a new, slower stream must begin its departure for the system 
to continue to satisfy the marginal stability condition. The velocity of this new 
stream will correspond to a new v lc0 . In the limit, this becomes a continuous 
process with vi co decreasing at a slower rate than v max as the density increases 
at lower velocities. 

Finding v lco and v max as functions of time is a tedious numerical exercise. 
However, it turns out that analytical solutions may be obtained for the set of 
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discrete times at which a given vi co subsequently becomes the relevant v max . 
This idea will be discussed in greater detail in the following paragraphs. 

The time, 


to 


1 

Vico 


( 3 . 42 ) 


at which a v, co flux stream becomes a v max stream is just equal to the time it 
takes the last electron in the vi co flux stream to travel from the other end of the 
kernel to the moving front and escape. In other words, the v lco stream becomes, 
momentarily, just before the last electrons pass through the front a v max itself. 
Using this key idea we can calculate a progression of i>j co ’s as a function of time in 
the following manner. Knowing v max initial, we first calculate the kernel length, 
L{t), at which time the last electron in the escaping flux stream reaches the front. 
This gives us a time, which the next V( co begins its departure from the 

kernel. Calculating the time of transit of the last electron in this flux stream 
provides the time at which this v lco becomes a v max itself. Thus, using this step 
by step approach we may calculate a sequence of vi co 's exactly. 

To calculate the length of the kernel when v lco {t = <„_i) becomes v max we 
write 


L 


n 


w 


kLnT n 


rrii 


’t n + L 0 



( 3 . 43 ) 


where t n is the specific time at which the last electron in the flux stream departs 
and the time when the particle reaches the moving front is 


Substituting, we obtain 


tn = tn - 1 + — . 


Vlco 


( 3 . 44 ) 
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and 



(3.45) 


. s „ fkLJ^ , , *> [^To L " , i 

L « ! = 3 \j~^r + 3 V m, v,„ +I, “ 

This takes the form (x = L*), 


x 3 + ax 2 + 6 = 0 


where 


= _ 3 fkL 7 k 1 
V m i v lco 


and 


jkLoTo j | 

6 = — 3a/ t n -i ^o- 

V m « 

T ak ing the only real root, we find that 


_i 

3 


(3.46) 


(3-47) 


(3.48) 


(3.49) 


x a a 2 / a 3 6 / a 3 6 6 2 \ 

^•=-3 + IlT-i + V¥ + 7j 


+,&?)*. (3.50) 

+ v 27 2 V 27 4 y 

Now, this value can be used in equation (3.44) to calculate the time at which the 
electron reaches the barrier. 

Using the marginal flux condition, the new v lco can be calculated from 


/*V m0 x 

J = a 2 / v 3 exp(-b 2 v 2 )dv 

Jv tco 


(3.51) 
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(3.52) 

(3.53) 


giving 


/ = 


3 

2 meV 


2 

e 


C s 


(3.54) 



(3.55) 


This can be solved to give a new vi co . As stated before, by setting each v/ co 
as a v max for the next step, the problem may be solved exactly at these times 
for any number of iterations. Appendix B includes the code which accomplishes 
these calculations. 


The initial condition flux rapidly departs. In other words, vi co decreases 
rapidly while v x , the edge of thermal bulk distribution in phase space, is increas- 
ing rapidly (Figure 3.8). At time — 1.4 seconds the two curves intersect (Figure 
3.9) implying that by this time the high-energy, initial condition electrons have 
escaped to the tail and that the kernel is completely thermalized. 

Thus, the MacDonald et al. (1957) analysis, while informative, is only 
needed to deal with the first second or so of the event. The front is main- 
tained thereafter by new high-energy electrons created by the pumping action of 
the collisional distribution. Hence v lco S* 2.6 v e (t) may be used to describe the 
period of time / > 1.4 seconds. 
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Figure 3.8z Hybrid Distribution 

This shows an example of the a distribution constructed from a Maxwellian 
kernel, a Maxwellian tail from a hotter distribution, and a spline linking function. 
It represents the distribution in the Model T kernel early in the event. Region I 
bounded at the higher velocities by t* is a cooling Maxwellian, region III is the 
collisionless tail remaining from the initial distribution, and region II, bounded 
at the higher end by t*, is a mixture of both. All flux streams above v lco axe 
escaping and all flux streams above v max have already departed. 
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Figure 3.9: Model T Velocity Comparison 

This shows the evolution of V[ co , the lower boundary of the flux escaping through 
the ion-acousic fronts, and t>i, the upper boundary of the Maxwellian bulk dis- 
tribution as a function of time. Clearly, the overlap occurring between 1 and 2 
seconds into the event suggests that at this point the kernel electron distribution 
has become entirely Maxwellian. 
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We now have reached the point in the development such that the electron 
distribution is known as a function of time for the entire event. Now using 
this distribution, both the thermal and nonthermal hard x-ray emission may be 
calculated for Model T. We use the Kramer cross section as a good approximation 
to the more rigorous Bethe-Heitler cross section (Li and Emslie, 1990; Emslie and 

Machado, 1987), 


where 


(3.56) 


k = 7.9 x 10 _25 Z 2 


(3.57) 


and 


Z 1 = 1.4, (3-58) 

the mean square atomic number for solar abundances (Allen, 1973; Emslie et al., 

1986a; and Tandberg-Hanssen and Emslie, 1988). 

Now, from Tandberg-Hanssen and Emslie (1988), the Maxwellian may be 

written as (electrons cm 3 erg ) 

f(E) = . 2ne — E*e~& (3.59) 

and the thermal emission at the Earth in photons cm keV s may be found 
from 


/,»(«) = jT fe(,E)v(E)c, r(e, E)dE, (3.60) 

where R is the mean distance from the Earth to the Sun. 
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This integral may be used to calculate the thermal bremsstrahlung from the 
bulk Maxwellian and the cubic spline regimes by integration from e to ej = \mv\ 
using a Maxwellian f(E) summed with an integration from Vi to v 2 using a cubic 
spline f(E). In a similar fashion, after 1 second, the thermal emission from the 
kernel may be found using a Maxwellian truncated at 2.6 v t . 

To illustrate this further, the intensity at the Earth due to a truncated 
Maxwellian for e < a is (photons cm -2 keV~ 1 s — 1 ), 


Ith(e) = (9.55 x 10" 13 ) 


rrlA 


e — e & 


(3.61) 


' eT 1 1 2 R 2 

where a = 3.38 m e i; 2 , / is the pixel length, R is the mean distance from the 
Earth to the Sim, and A is the flare cross sectional area. For e > a, we have 


T th = 0. (3.62) 

For the nonthermal emission, we follow the approach used in Chapter II 
with the difference that the injection point (N = 0) will now correspond to the 

time- varying edge of the kernel. For a Maxwellian tail (Brown and Emslie, 1988), 
we have 


F(E 0 ) = 


2t n e E 0 

(irm e )2 (A :T)§ 


kT 


(3.63) 


in electrons cm 2 s -1 erg -1 . 

Using equation (3.63), the cumulative intensity may now be found numerically 
from (Appendix B) 


J( e i N) = — f N f 
e Jn '= o Jb 


F(E 0 ) 


(E 2 - 2KN) 1 / 2 
where b = max[E c , (e 2 + 2KN)*] and E c = 3.38 m e u 2 


dE 0 dN' 


(3.64) 
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The nonthermal emission in photons cm 


-2 s -i keV -1 for the i’th pixel may 


be found from 


Ii(e) = J(e,Ni) — 

Figures 3.10-3.13 show the results of this analysis for photon energies equal 
to 10, 40, 70, and 100 keV, respectively. At 10 keV (Figure 3.10), the apex is 
extremely bright with weaker emission from the footpoints. The effect of the 
growing, cooling kernel can be seen in the emission as well. Some nonthermal 
emission is found between the kernel and the footpoints but is dominated by the 

kernel and footpoint emission. 

As one goes higher in energy, the emission from the apex decreases much 
faster than that in the footpoints until at 40 keV (Figure 3.11) the apex and 
footpoints are emitting roughly the same flux at t = 0. At 70 keV (Figure 3.12) 
and 100 keV (Figure 3.13), the apex is still visible but the emission from all 
sources is greatly reduced from that at lower energies. Indeed, at higher energies 
the Model T mimics the emission profiles of the Model N-T discussed in Chapter 

II. 


C. SUMMARY 

This chapter developed the Model T quantitatively to provide emission pro- 
files as functions of time, energy, and loop position. By examining the results 
of this model and comparing them with the Model N-T from Chapter II we im- 
mediately ascertained that differences in the emission profiles of the two models 
were indeed present at lower energies. At higher energies, the Model T profiles 
began to mimic those of the Model N-T. 
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Figure 3.10: Model T Emission Profile At 10 keV 

This shows the emission for 21 pixels composing the loop. Early in time the 
bright central pixel is dominant. As the kernel expands and cools, the central 
pixel emission decreases while that of adjoining pixels increase. In general, the 
overall signature is much different from that of the Model N-T at the same energy. 



Figure 3.11: Model T Emission Profile At 40 keV 

This shows the emission for 21 pixels composing the loop. In this profile early 
in time, the central pixel is as bright as the footpoint pixel. Due to the rapid 
cooling of the kernel, the emission from both pixels decreases rapidly. In general, 
after one second, this signature begins to look much like that of the Model N-T 

at the same energy. 
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Figure 3.12: Model T Emission Profile At 70 keV 


This shows the emission for 21 pixels composing the loop. In this profile early in 
time, the central pixel is almost as bright as the footpoint pixels. This is expected 
as the highest energy electrons escape from the kernel to produce emission in the 
footpoints. Except for the short burst of emission from the apex at t = 0, this 
signature, in general, resembles that of the Model N-T at the same energy. 
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Figure 3.13: Model T Emission Profile At 100 keV 

This shows the emission for 21 pixels composing the loop. Very early in time, 
the a three-point signature is present; however, in a short time a situation is 
reached where relatively few high-energy electrons are available in the kernel for 
producing high-energy photons since most have escaped to the footpomts. In 
general, after t = 0, this signature will resemble that of the Model N-T at the 


same energy. 


Chapters IV, V, and VI will discuss the two telescopes which we shall use to 
view these profiles and examine the capabilities of these instruments. Chapter 
VII will meld together these chapters and the results from Chapter II and III 
to determine if these profiles can indeed be imaged with a Fourier telescope to 
distinguish between the two models. 
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IV. THE FOURIER TELESCOPE 


A. INTRODUCTION 

Chapters II and III discussed the emission from Model N-T and Model T as 
functions of time, energy, and position. A detailed comparison will be discussed m 
Chapter VII; however, it is sufficient here to point out that discernable signatures 

were discovered to exist between the two models. 

Model N-T was characterized by a footpoint structure which subsequently 
evolved into a filled loop. By contrast, from Chapter III, Model T initially 
exhibited a bright central region in addition to the footpoint emission from the 
escaping high-energy electrons. This central region spread downward along the 
loop as a consequence of the diffusion of heat into the surrounding plasma. These 
spatial emission signatures are in principle sufficient to allow a discriminating test 

between the two models through imaging. 

Now, the question becomes whether a Fourier telescope can view those pro- 
files to provide the observer with a discriminating set of images. This chapter will 
discuss the history and introduce the basic concepts which form the foundation 
for Fourier telescopes while Chapters V and VI will discuss the spatial modula- 
tion collimator and the rotating modulation collimator, respectively, in greater 
detail. Chapter VII will then discuss the capabilities of the Fourier telescope to 

image the signatures derived in Chapters II and III. 

Hard x-rays (10-100 keV) cannot be imaged by conventional optics such 
as lenses or grazing incidence mirrors. Only recently has the technology been 
developed to allow these types of sources to be imaged with reasonable resolu- 
tion. Surprisingly, this technology is an extension of the techniques used in radio 

astronomy. 

Essentially, the Fourier telescope concept involves sampling selected Fourier 


79 


components from a wavefront emitted by an extended source on the Sun’s surface. 
Each component may be visualized as a single point on a common complex 
surface. By measuring a number of discrete components over a sufficiently large 
spatial frequency spectrum, this Fourier surface may be approximated. A Fourier 
transform of this surface function yields an approximate or a dirty image. The 
heritage of radio astronomy comes to the rescue in that several algorit hms , which 
have been proven over the years, have been developed to clean the dirty image 
to produce a more meaningful result. 

Since obtaining final images requires a tremendous amount of numerical 
processing, end-to-end numerical simulations are vital to both the design and 
the data analysis phases of a mission. In the design stage, the simulations can 
help keep track of performance gains or losses versus engineering trade-offs, and, 
m the data analysis stage, the observer needs to know whether the image he has 
obtained was created with enough photons to insure good imaging. 

In principle, Fourier telescopes can be divided in two general classes: those 
that rotate (rotation modulation collimators) and those that do not, spatial 
modulation collimators (Campbell et al., 1991a; Hurford,1977). 

Figure 4.1 illustrates a basic telescope geometry with two grid planes fol- 
lowed by a detector assembly. Grid pairs can be constructed to modulate the 
incoming wavefront over the detector allowing a particular component to be 
measured (Figure 4.2). Typical geometries from proposed telescopes were used 
to develop the telescope models described in Chapters V and VI. Instrument 
imperfections such as twist, bending, uneven thermal expansion, etc. can als o 
be addressed individually and in combination using simulations. 

As mentioned earlier, conventional imaging techniques use reflection or re- 
fraction of the incident radiation (e.g., Hecht and Zajac, 1976). Since these 
processes do not work well at wavelengths < 2 A, non-focusing Fourier colli- 
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Figure 4.1: Basic Fourier Telescope Geometry 

This shows a conceptual Fourier telescope. Each plane is composed of several 
grids with slit spacings corresponding to the spatial frequencies to be measured. 
The bottom plane’s grids are directly related and aligned to those in the top 
plane so that each Fourier component is provided by a gnd pair illuminating a 

detector. 
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mation techniques must be used for image formation. The two-grid modulation 
collimator first introduced by Oda (1965) provided a single Fourier component. 

It had the advantage of simultaneous measurements of signal and background; 
and, it also had the ability to locate a source. 

Extending the concept (Bradt et al., 1968) to multiple grids or rotational 
s cannin g permitted the formation of true images, either continuously spread 
across the image plane or subdivided into discrete picture elements by subcolli- 
mators that feed independent detector elements (van Beek, 1975). A practical 
limit on the angular resolution of these multiple-grid collimators was the need 
for mechanical rigidity of the support structure. The structure had to maintain 
the relative positions of the grids to a fraction of the individual slit width which 
in practical systems (van Beek, 1975) may be as small as 50 /x. 

Makishima et al. (1978) developed a two-grid modulation collimator that 
measured the corresponding Fourier component of the angular distribution of 
the source. Thus, a single subcollimator determined the identical parameter 
measured by a two-element interferometer as used in radio astronomy. As m 
aperture synthesis (Fomalant and Wright, 1974), a number of such measurements 
can be combined to produce an image. 

B. Theory 

Consider a point source located some distance from the observer. The emit- 
ted radiation may be envisioned as an electromagnetic spherical wave emanating 
from the source into space. When this wave reaches the instrument in low-Eaxth 
orbit it is essentially planar over the extent of the telescope (typical diameter of 
1 meter). Mathematically, this wave could be described by taking the Fourier 
transform of the brightness distribution of the point source (Steward, 1987). In 
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two dimensions, the Fourier transform for the brightness distribution B(x,y) 
may be written generally as 

/ oo 

B(x,y) e 2 ^ ux+v ^dxdy (4.1) 

■OO 

where u and v are spatial frequencies. Taking the inverse transform, 

B(x',y')= f f(u,v)e- 2iri ( ux ' +v y''>dudv (4.2) 

J —oo 

provides an image. 

Most optics texts discuss extensively the use of the Fourier tr ansf orm in a 
related operation in which the E field vector associated with an electromagnetic 
wave emanating from the point source and measured a great distance from that 
source will be transformed. This transformation may be described mathemati- 
cally by a Fourier transform and the intensity in the image plane may be obtained 
by taking the modulus squared of the complex vector. This particular approach 
is especially useful in treating diffraction and many other wave-related effects. 
There is an associated approach in which the brightness distribution from the 
source as measured a great distance away may also be represented by a Fourier 
transform. Taking the inverse transform in principle provides an image directly. 

In the laboratory an analogous experiment could be performed in which 
a source was placed at the left focal point of a lens. Another lens would be 
placed behind the first and at the right focal point of this lens an image could be 
observed. However, if the intensity distribution is measured between the lenses, 
one finds that the radiation has structure in this region as well from which a 
legitimate image may be obtained. 

As mentioned above, there is an interesting but not altogether unexpected 
parallel with radio astronomy which also constructs its images from measured 
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Fourier components. The two-element interferometer provides high resolution 
by correlating the signals of the two antennas (Figure 4.3). The correlation is 
normally achieved by the multiplication or addition of the signals, which produces 
a spatial modulation of the primary beam of the antennas with interference 
fringes (Figure 4.4). In this way fine structure is introduced into the primary 
beam to increase the resolution. 

The response of the system to a point source of monochromatic radiation 
of frequency w or wavelength A can be described by a voltage E proportional to 
the sum of the electric fields generated at the feed of each telescope at slightly 
different times. The time difference is called the geometric delay and is denoted 
by t. The voltages at the multiplier input are 


V\ oc E cos(u>t) 

(4.3) 

y 2 oc E cos[uj(t — t)] 

(4.4) 


and 

y 2 = E cos{ut - - y-cosd ) ( 4 - 5 ) 

where B is the separation of the two antennas, A is the wavelength of the ra- 
diation, and 0 is the angle between the point source and the line joining the 
two antennas. The expression (B/ A) cosO , the interference term, gives the phase 
path-length difference of the radiation travel along the two possible paths. The 
diurnal motion of the Earth causes 6 to vary with time. The output, R(t), of the 
multiplier, after a high-frequency term is rejected by a low pass filter, is 

R(t) oc S cos^y~ cos6(t )] ( 4 - 6 ) 


85 


Diurnal 



R It) 


Figure 4.3: A Basic Radiotelescope 

The illustrates the concept of a radiotelescope. The principle here is that the 
time difference in the arrival of the wavefront leads to an interference pattern 
which is directly related to a Fourier component. 
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Figure 4.4: A Typical Interferometer Response 

The envelope of the response is sinusoidal in nature and directly related to a 
Fourier component in that one may obtain a phase shift and an amplitude which 
are sufficient to define a phasor in the complex plane. The real and imaginary 
components of the visibility function are then easily obtained as components of 

the phasor. 
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This is the basic equation of interferometry. The flux density, or power S, of the 
source has replaced E 2 . The fringe spacing is given by the angle which produces 
a change of one wavelength in the path-length difference. 

The phase path- length difference (y) cos 6 can be more generally written 
as B • s where B , the physical spacing, is equal to the element separation in 
wavelengths and its direction is that of the line joining the elements. The baseline 
direction will be defined as a vector from telescope 2 toward telescope 1. The 
preceding equation then becomes 

R(t) a S cos[ 2 ttB • $(<)] (4.7) 

The response of a two-element interferometer to an extended source can 
be obtained by considering the source to be a collection of point images and 
summing their individual responses. Let s be the vector to a convenient position 
near the source. This point is denoted as the phase center. Any other point 
can be denoted by s + a. If /(<?) describes the brightness distribution (angular 
distribution of power), then the response to the extended source is 

-R(<) = J I{&) cos[ 2 irB -(s(t) + 0)\d&. (4.8) 

A useful quantity in interferometry is the projected spacing b of the physical 
baseline B as viewed from a radio source. The projected spacing b is given by 

b = s x (B x s) (4.9) 

= B-(s-B)s. (4.10) 

Generally, the projected spacing is resolved into components along directions 
to the east and north, which are commonly denoted u and u, respectively. 
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Since the angular size of the region observed is limited by the extent of the 
antenna response (typically less than 1°), the phase term can be expanded to 
first order for sufficient accuracy, 


B-(s + d) = B-s + B-a 


(4.11) 


m 


« B ■ s + bd ( 4 - 12 ) 

Since a is nearly perpendicular to s, only the projected spacing b is retained 
the second term of the cosine. The response becomes 


/ +°° _ 

1(d) cos[2irB ■ s + 2-irb • a] 

-OO 

This may be expanded into 


da 


(4.13) 


R(t) 


= cos(2irB s) [ 1(a) cos(2ivb-a) da (4.14) 

J — OO 


— sin(2-KB ■ s) / 1(a) sin(2nb ■ a)da (4.15) 

J — OO 

However, it is much easier to work with the more compact complex form 


/ +OO 

KS^^-’dS). (4.16) 

-OO 

Now, define the visibility function by the following expression, 


V = t I(a)e 2 * its da. (4-17) 

J — OO 

This integral is a complex number. The amplitude of V is proportional to 
the amplitude of the fringe pattern and the argument of V equals the phase 
shift in the fringe pattern from that of the response to a point source at the 
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phase center. Clearly, the visibility function is the Fourier transform of the 
brightness distribution. Since, in radio astronomy, the time variability of the 
source is normally less than that due to the diurnal motion of the source, the 
instantaneous response can be rewritten as (Hurford and Hudson, 1980) 

/ +oo r+oo 

/ I(x,y) e 2 ^ x+v ^dx dy} (4.18) 

“OO J — oo 

where K is independent of the source. The integral has been explicitly written 
in terms of eastward and northward displacements (x, y) on the sky, and (u,v) 
are the corresponding components of b(t). Note that the response depends only 
on the spatial frequency content in the source that corresponds to u and v . 

This can be further simplified by rewriting the integral in terms of (0, <f>), 
orthogonal coordinates in the sky parallel and perpendicular to b , 


R = i?e{ K j 

r+0max 

f F(6) e 2nibe dd) 

&m i n 

(4.19) 

where 

/•4'^max 


m = 

/ WA)d4> 

(4.20) 

is the projection of the source brightness distribution perpendicular to b. 

To characterize the modulation collimator response (i.e. Moire fringes) in 

x-ray astronomy, the transmission pattern may be represented by 
wave A(0) that repeats with angular period P: 

a triangular 

m 

_ , m 

p ’ 

(4.21) 

for 



-P 

^ p 


2 

< 6 < — 
~ ~ 2 

(4.22) 


and 


A (6±P) = A (d). ( 4 - 23 ) 

A scanning modulation collimator sweeps out the angular range P m a 
finite time, and the counting rate as a function of time must be converted to 
the angular profile from an aspect solution that describes the collimator motion. 
The angular response can be modeled with a function N{8 ) noting that N(8)A8 
is the probablility of AN counts in the angular range A 6: 

/ + Ipmai / A .v 

[0.5 F(ip) A (6 — VO + B] At dip, (4.24) 

'ifimin 

where F(ip) is the source brightness contribution (photons cm -2 s -1 radian 1 ) in- 
tegrated parallel to the collimator aperture (a 50% maximum reponse is assumed 
for slit widths and boundary material of equal spacing); B is the background rate 
in counts cm" 2 s" 1 , A the detector area (cm 2 ) and At the integration time (s). 
xp is the brightness distribution parameter across the grid. 

The complex visibility function as measured by a given collimator may be 

written as 

v = 1 j % N(«) e i * e id. (4-25) 

Substituting, the following is obtained 

V = — f^ maz f 2 A (0 -tp) + 2 B] At e' 2 ? 0 dd dip. (4.26) 

P J-Tpmin >'-■§ 

Since A (8 - ip) is periodic, this becomes 

v = AA( f + *— FWe i^ P iff. (4.27) 

P J — lpmin 
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If we now approximate the periodic function with a sinusoid with the same am- 
plitude and frequency and substitute for A(0') we obtain 



where c\ is a constant. 

This has the same general form as the response for the radio interferome- 
ter examined previously (equation 4.17). However, we note that the triangular 
response used in this analysis offers in principle higher Fourier harmonics which 
may potentially be used to provide higher resolution information about the source 
(Hurford and Hudson, 1980). 

The implication of this correspondence for x-ray astronomy is that it may 
take advantage of many of the techniques and methods which have been rigor- 
ously developed over the years for radio astronomy. This identity is also interest- 
ing in that for radio waves the characteristic periodic variation of response with 
angle comes from wave interference, while for x-rays it comes from geometrical 
shadowing. 

Along these lines, a Fourier transform telescope consists of a set of individual 
modulation collimators (sub collimators) that feed discrete detectors. A subcol- 
limator is equivalent to one baseline or a simple two-element interferometer in 
radio astronomy. A measurement of the amplitude 


V = \/[f2e(y)] 2 + [Im(V)]2 (4.29) 


and phase 


a 



(4.30) 
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determines all of the necessary information regarding the Fourier component 
measured by a given sub collimator. 

Re and Im are the real and imaginary representation of the components of 
the complex visibility function in the complex plane at a single point in the (u,v) 
plane. Since the visibility function V(u,v) in the observational plane is directly 
related to the brightness function in the object plane via a Founer transform, 
an image is obtained from the inverse Fourier transformation of the visibility 

V(u,v). 

The corresponding (u, v) coordinates are just 

U = lease ( 4 - 31 ) 

and 

v = ^ sin 6. (4-32) 

where P is the angular period of the grids. For slits and solid slats of the same 
width, P = (2 s/d), where s is the slit/slat width and d is the grid separation. 
The collimator period may be related to (u,u) by 

p = = — (4.33) 

+ v 2 

Finally, it should be noted that Fourier telescopes are photon counting in- 
struments. In other words, the telescope takes a number of photons and through 
numerical reconstruction provides an integrated image. It leaves to the observer 
the task of organizing the photons into energy and time bins. Usually, the ob- 
server accomplishes this with an innovative detector design. This is a powerful 
feature of this type of instrument since, in the case of insufficient flux, it allows 
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one simply to combine bins, both in energy and/or in time, to provide enough 
photons for an image. 


C. THE ASTRONOMICAL IMAGING PROCESSING SYSTEM (AIPS) 

AIPS, developed and maintained by the National Radio Astronomy Ob- 
servatory (NRAO), is widely accepted today in the radio astronomy commu- 
nity (Greisen, 1986). It is a highly sophisticated software system requiring over 
180,000 blocks of VAX disk space to run. At the core of the multitude of software 
routines are those which allow the user to bring Flexible Image Transport Sys- 
tem (FITS) formatted (u,t>) data into AIPS (Greisen and Harten, 1981; Wells, 
Greisen, and Harten, 1981); fast Fourier transform the data to form a dirty map; 
clean the dirty map in a number of independent ways to form an image; and dis- 
play the image in a number of different ways. Peripheral software allows the user 
to enhance an image by rotating it, zooming in or out, or suppressing background 
noise. Figure 4.5 shows the AIPS user manual cover. It actually does contain 
recipes such as banana nut bread in addition to being an excellent reference. 

AIPS uses the equatorial system of coordinates. In this system, the Earth’s 
equator is the plane of reference. The poles are at the intersection of the Earth’s 
axis with the celestial sphere, an imaginary surface at a large distance with the 
Earth at its center. Thus, the ordinate of AIPS output images will be given as the 
decimation m units of degrees, arc minutes, and arc seconds while the abscissa 
will be given as the right ascension in units of hours, minutes, and seconds. 

One major advantage in using AIPS is that it allows a dirty image to be 
cleaned in several different independent ways. Two of the most widely used 
are the CLEAN algorithm and the Maximum Entropy Method (MEM). In the 
CLEAN method the respective phases and amplitudes are combined with the 
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Figure 4.5t AIPS Cookbook Cover 

(NRAO, Virginia, 1986) 
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measured Fourier components to produce an uncleaned image. This preliminary 
image generally has large sidelobes, particularly if there are one or more dominant 
point sources of emission in the field of view. 

The CLEAN algorithm attempts to reduce these sidelobes by iteratively 
deconvolving the image with the expected response of an ideal point source in- 
cluding sidelobes. This method has been shown to be effective for a wide variety 
of source configurations, although the method works better for collections of 
point sources than for diffuse emission regions. The advantages of the CLEAN 
method are that it is relatively fast and theoretically simple to apply. Some dis- 
advantages are that the method can produce images with areas of negative pixel 
brightness, that statistical errors are somewhat difficult to incorporate, partic- 
ularly when the errors vary from component-to-component, and that significant 

sidelobes can still be present in the image, particularly if the statistical errors 
are large. 

The Maximum Entropy Method attempts to produce the smoothest (highest 
entropy) image compatible with the data. Advantages of the method are that 
all image pixels have positive brightness and the statistical errors of individual 
Fourier components are easily incorporated into the analysis. While MEM has a 
bias toward diffuse sources of emission, in practice the method works well even 

for point sources, producing images with lower sidelobes and reduced artifacts 
than the CLEAN algorithm. 

One of the strengths of the numerical model is that a dirty image set can be 
cleaned conveniently using both techniques and then the images compared. By 
using known inputs, the techniques can be compared in different scenarios and the 
best one chosen. In addition, artifacts which might be introduced by assumptions 
in the techniques will become apparent and will be useful information when real 
data are obtained and processed through the model. AIPS, obtained from NRAO, 
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has been adapted to process x-ray data. 

As part of the verification process, (u, v) data were obtained from Palmer 
and Prince (1987). In this experiment, a Fourier telescope was developed in 
the laboratory as a pathfinder for future development. Grids were fabricated 
from soldering wire, and small bits of radioactive cobalt were used as hard x- 
ray sources. Since a steady state source was used, the telescope could view as 
long as was necessary to obtain enough photons for an image. Also, the source 
distribution could be rotating rather than the telescope. The experiment was a 
success and paved the way for more advanced developments. 

Figure 4.6 shows the Palmer and Prince (1987) results using three dimen- 
sional contour plotting, gray scale graphics, and contour plot images using both 
CLEAN and MEM. While the author noted that maximum entropy seemed to 
produce a higher quality image, no such difference was found here in his data 
using AIPS. This may be attributed to the fact that the AIPS algorithms are 
highly sophisticated and under scrutiny from many users and are constantly be- 
ing improved. However, for more complicated images such as the filled loop 
(Chapters V and VI), MEM did provide better images than CLEAN. Figures 
4.7, 4.8, and 4.9 show the AIPS output using the Palmer data in contour, three 
dimensional, and grey scale. The results are reproduced exactly with the AIPS 

CLEAN algorithm. 

In order to work with these data, a complicated linking program was required 
to convert the FORTRAN data into a FITS file acceptable to AIPS. The data 
were then brought into AIPS and fast Fourier transformed to obtain a dirty 
map. This map was then cleaned using several different algorithms and displayed 
in several different graphic formats. AIPS duplicated the results found in the 
literature exactly, thus greatly increasing our confidence that the software had 
been successfully implemented. 
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FIGURE 5. GRAY SCALE PLOT OF MAXIMUM FIGURES. THREE-OIMENSIONAL CONTOUR PLOT 

ENTROPY IMAGE SHOWN IN FIGURE 48. OF THE MAXIMUM ENTROPY IMAGE SHOWN IN 

FIGURE 48 AND FIGURE 5. 



FIGURE 4A. CONTOUR PLOT OF THE IMAGE OF FIGURE 4a CONTOUR PLOT OF THE IMAGE OF 

THREE 57 Co SOURCES PROOUCED USING THE THREE 17 Co SOURCES PRODUCED USING THE 

CLEAN ALOGRITHM. CONTOURS AND SIZE OF MAXIMUM ENTROPY METHOD. THE DATA SET 

FIELD OF VIEW ARE THE SAME AS IN FIGURE 3A. FOR THE IMAGE WAS IDENTICAL TO THAT USED 

TO PRODUCE FIGURE 4A. 


Figure 4.6: Palmer Results 


Palmer built a laboratory version of a spatial modulation collimator using sol- 
dering wire for grids and three cobalt sources. His telescope was successful in 
producing in the laboratory hard x-ray images of the sources. 
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Figure 4.7: AIPS Contour Of Palmer Results 

Using data from an actual laboratory spatial modulation collimator, AIPS du- 
plicated the Palmer’s results. The contours in this image representation are 
percentages (10 % increments) of the peak value in the field of view. 
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Figure 4.8: A IPS 3-D Of Palmer Results 


It is useful to have both three-dimensional displays and contour displays with 
which to analyze data. AIPS provides both. Here, AIPS has processed data from 
a laboratory telescope to duplicate those results. 
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Figure 4.9: AIPS Grey Of Palmer Results 

This is another example of the versatility found in AIPS. Here, using data 
from a basic laboratory telescope, a duplicate grey scale image was obtained. 
The contour in the center is a representation of the point spread function of the 

telescope. 
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D. POTENTIAL ENGINEERING CONSTRAINTS 


Numerical models typically precede the actual hardware development of a 
complex system and generally if the system does not perform well numerically 
then its development is terminated at that point. Hence, numerical models pro- 
vide relatively inexpensive proof of principles tested in the systems development 
process. 

In building the actual system the designer may be forced to accept compro- 
mises to the ideal due to cost and/or schedule constraints. Reducing the grid 
thickness is one example in which weight, cost, and manufacturing time might 
be reduced. The technical implications of such a change would be to allow high 
energy photons to penetrate one or both grids. These photons would essen- 
tially then constitute a source of random noise on the detector. Other sources 
of random noise would be the cosmic background and the detector itself (inter- 
nal noise). Fortunately, all of this results in a relatively uniform background of 
counts while the modulation provided by the telescope rides upon this uniform 

background level. This feature makes the Fourier telescope tolerant to 
random noise. 

Real detectors will offer other design challenges such as decreasing efficiency 
at higher energies, variable sensitivity as a function of time, and manufacturing 
difficulties for large sensitive areas. For our particular sytems operating from 10 
to 100 keV, efficiencies available from off-the-shelf detectors are close to 100%. 

The hardware development of the Fourier telescopes described in Chapters 
V and VI will have some risks associated with it and these risks will grow as 
one attempts to extend the technology to energies above 100 keV. For example, 
grid manufacturing still remains a high risk issue, especially for the finer ones 
providing 1-4 arc second spatial resolution. Figure 4.10 shows the throughput of 


102 


tungsten for different thicknesses. Clearly for 10-100 keV the grids can easily be 
made opaque using a modest thickness; however, as one moves to higher energies 
increasing grid thickness may be necessary. For the finer grids, this may present 

a formidable nianufactoring challenge. 

Given that these grids can be successfully built, they must then be mounted 
into a structure 5 meters long in such a way as to provide high precision alignment 
between grid pairs. The alignment error must be on the order of a few microns. 
The telescope must then be accurately pointed. 


Detector saturation is another high risk area. Around 10 ji, there is a danger 
of saturating the detector with too much flux and at energies above 100 keV 
detector efficiency may change drastically with energy. The Burst and Transient 
Source Experiment (BATSE) is presently in orbit using large, Nal(Tl) detectors 
(127 cm 2 ) in the 15 keV-110 MeV range for sensing hard x-rays and gamma rays. 
These compare very well with our 100 cm 2 detectors used in our two telescopes. 
From 10-100 keV, the BATSE detector efficiency is uniformly 90-100%. 

In addressing the saturation problem, Figure 4.11 illustrates the throughput 
characteristics of aluminum. Clearly, an aluminum shield of modest thickness 
covering the end of the telescope could be contructed to eliminate low-energy 
saturation problems. Indeed, by simply varying the thickness one could tune the 


shield to provide practically any desired response. 

Figure 4.12 illustrates for an actual bright flare the BATSE detector sat- 
uration and the SMC imaging thresholds. Clearly, for a thin aluminum shield 
0.04 centimeters thick we can image the flare from 10-100 keV without fear of 


saturation. 

Therefore, while our simulations will represent a best case for 
Fourier telescope performance, it will be a meaningful best case in 
that both grids and detectors can be reasonably constructed to provide 
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Figure 4.10: Tungsten Throughput 

For a very reasonable tungsten thickness of only 0.3 centimeters, the grids 
be made opaque to x-rays from 10-100 keV. 
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Figure 4.11: Aluminum Throughput 

This shows the throughput for different thicknesses of aluminum as a function of 
energy. Clearly, an aluminum shield of reasonable thickness can be constructed 
to eliminate flux below 10 keV and avoid detector saturation problems. 


105 



Figure 4.12: BATSE Detector And SMC Imaging Thresholds 

Using the Burst and Transient Source Experiment (BATSE) detector (currently 
operating in orbit) as a reference, and choosing a bright flare, clearly for a very 
reasonable thickness aluminum shield of 0.04 cm, detector saturation could be 
avoided completely while retaining a flux level well above the Fourier telescope 
imaging threshold. 
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close to ideal performance 


between 10-100 keV. 


E. SUMMARY 


Fourier telescopes can be designed in several different ways to optimize avail- 
able resources. The rotation modulation collimator approach, which will be taken 
in Chapter VI, will essentially consist of modulating the x-ray wavefront in time 
with two different grid pairs (i.e., modulation periods) corresponding to angular 
resolutions a 4 and 8 arc seconds. The temporal modulation will be accomplished 
by rotating the telescope, consisting in this case of 48 grid pairs representing the 
two spatial frequencies, about its line of sight axis. Each spatial frequency will 
be represented by 12 cosine and 12 sine grid pairs with a relative orientation 
between each of 15 degrees. This design will allow a relatively large number of 
Fourier components to be measured per grid pair and allow the use of a single de- 
tector rather than a spatially discriminating detector behind each. Furthermore, 
snapshot images of very bright sources may be obtained. 

Since the RMC must integrate over the time of rotation of the telescope, 
rotating the telescope faster may improve the temporal resolution. On the other 
hand, this will decrease the number of effective counts or photons reaching the 
detector for a given Fourier component. For solar flare processes which evolve 
over time scales of several minutes, this does not present a problem; however, 


most impulsive flare processes evolve over time scales of seconds. Also, for some 
basic RMC designs, the observer must have a priori knowledge of the rough 
location of the source in order to remove all ambiguities from the final image. 
Image information may also be lost as the source moves into the near vicinity 
of rotational axis of the telescope. All of these concerns raise the question as to 
how well flares may be resolved at the edge of the RMC performance envelope. 
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The importance of numerical modeling for exploring the validity of low photon 
count images as well as an optimization tool is evident. 

An alternative method which avoids the temporal resolution problem ex- 
perienced with the rotating collimator is to use spatial collimation. If a grid 
geometry is selected so that the bottom grid slits are a little narrower than those 
in the top grid, a grid pair can be designed to distribute the incoming beam in 
one period of a triangular response across the detector area. This spatial trian- 
gular response is equivalent to the temporal triangular response produced by the 
rotating collimator. Thus, a simple imaging collimator with one grid pair with a 
position-sensitive detector may measure one Fourier component of the source 
angular distribution. Unfortunately, this Moire pattern now requires a spatially 
sensitive detector array to measure intensity across the pattern. Thus, the trade- 
off m going this way is more inherent engineering complexity as compared to the 
RMC. 

Chapters V and VI will discuss the spatial modulation collimator and the 
rotating modulation collimator in greater detail and Chapter VII will discuss a 

comparison between images of synthetic flares constructed using Models N-T and 
T. 
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V. THE SPATIAL MODULATION COLLIMATOR 
A. INTRODUCTION 

In Chapters II and III, significant differences in emission profiles as a function 
of time, energy, and loop position were found to exist between Model T and 
Model N-T. This chapter will describe the construction of a spatial modulation 
collimator which will later be used to view these profiles to determine if the 
differences can be imaged. Chapter VII will discuss the differences between the 
images in greater detail. In addition, some basic Fourier telescope degradations 
and limitations will be discussed in this chapter and in Chapter VI. 

A basic telescope design using 48 grids, but only two spatial frequencies, was 
employed for this design. Design parameters were chosen to be representative 
of instruments being flown at the current time and proposed to be flown in the 
near term. The Fortran code used to simulate the telescope may be found in 
Appendix C. Figure 5.1 illustrates the (u,v) plane coverage for the telescope. 
Note the limited coverage afforded by this design; yet, as we will see later 
Chapter VII, the telescope shall perform reasonably well against the physical 
profiles from Chapters II and III- 

This particular telescope uses a grid geometry in which the bottom grid slits 
are a little narrower than those in the top grid such that the incoming beam is 
distributed in one period of a triangular response across the detector area. If 
the grids are copied onto acetate and then superimposed, one obtains a smgle 
period Moire fringe pattern equivalent to the detector pattern produced by a 
single point source. Moving the grids apart in a direction perpendicular to the 
slits is equivalent to moving a point source away from the telescope line of sight. 
Figure 5.2 and Figure 5.3 illustrate the principle. By moving Figure 5.2 over 
Figure 5.3, one can generate Moire fringe patterns similar to those which would 
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Figure 5.1: SMC Fourier Plane Coverage 

The apparent limited coverage of the (u,v) plane is illustrated here. This gives 
rise to noise internal to the processing which manifests itself as artifacts in the 
final image. Fortunately, techniques such as the maximum entropy method can 
be used to reduce this noise. 
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be seen by an x-ray detector. 

Since we axe measuring broad fringes, only one-dimensional, low resolution 
detectors are necessary. Figure 5.4 shows the effect of a point source located 20 

arc seconds from the telescope pointing axis. 

Each photon as its strikes the detector along the local * arris may be viewed 
as a phasor in the complex visibility plane with amplitude one and a phase angle 
given by its x position. Since x varies from 0 to 2rr, the phase may be read directly 
from the data for that particular photon. A resultant may then be calculated 
and the real and imaginary parts of the Fourier component may then be found 

from 


and 


n 

= cos ^ i 
i= 1 


(5.1) 


n 

Jm = ^2 
»=i 


(5.2) 


where n is the number of photons reaching the detector. The associated («,») 
coordinates may then be found as described earlier in Chapter IV. 

The generation of this Moire pattern now requires a one-dimensional, spa- 
tially sensitive detector array measuring intensity across the pattern. Palmer 
and Prince (1987) determined that the pixel sire of the detectors need only to 
be of modest sire. This is because we are taking advantage of the fringes pro- 
duced by the interference of the grid pair with the incoming wavefront. As found 


elsewhere in optics and radio astronomy, these fringes are extremely sensitive 
to the location of the source. A point source at any angle will produce a peak 
counting rate somewhere on the detector; the location and height of this peak 
corresponds to the phase, *, and the amplitude, A. Thus, one grid pair with 
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Figure 5.3; SMC Lower Grid 


This grid contains (n + 1) slits and slats which in conjunction with the (n) top 
grid produces a 1 cycle beat pattern on the detector (Crannell, 1988). 
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PRECEDING BLANK NOT FILMED 


Figure 5.4: SMC Moire Fringes For Point Source At 20” 


Note the fringe pattern is somewhat similar to a radio telescope response with 
a fine modulation pattern overlying a more coarse pattern. The coarse pattern 
contains the information about the source. The location of the peak of the pattern 
from the center of the grid and the amplitude of the pattern are s uffi cient to 
provide both the real and imaginary components of the visibility function. The 
slit width, grid pair separation, and grid orientation provide the (u, v) coordinates 
for these components. 


a position-sensitive detector can measure one Fourier component of the source 

brightness distribution (Campbell et al., 1991b). 

The combination of data from many such subcoffimators, each with its 
position-sensitive readout, permits the synthesis of an image (Crannell, 1988). 
Different slit orientations and spacings generate the different Fourier compo- 
nents needed. In essence the Fourier-component measurement consists of a three- 
parameter fit (A, *, n) to each subcollimator output, where n is the mean count 

level. 

The geometry used for our SMC includes two slit widths of 0.0275 and 0.0125 
cm with a 500 cm separation between grids. This allows the measurement of two 
spatial frequencies over a set of equally spaced orientations in the (u,») plane. 
Using this telescope geometry, an end-toend Monte Carlo model (Appendix C) 
was developed using the spatial modulation collimator approach to investigate 
telescope performance under various flux conditions. Individual detector sire was 
selected to be 100 cm 2 for a total sensitive area of 4800 cm 2 . Table 5.1 compares 
our SMC to the Solar-A (YOHKOH) and HESP instruments. 

B. THEORY 

In principle a single high intensity burst of photons from an extended source 
will provide exactly the same image as the same source emitting a single photon 
sequentially from random points within the source provided the population of 
photons in both cases is large enough for imaging. This principle is the key 
to allowing one to simulate optical systems numerically. Of course, in practice, 
one is limited by the performance characteristics of real detectors (e.g., detector 
saturation). However, as discussed in Chapter IV for the energy range 10-100 
keV, we may be confident that our results are representative of the best case 

performance for the telescope. 
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SMC 

YOHKOH 

HESP 

grid material 

tungsten 

tungsten 

tungsten 

grid pairs 

48 

64 

12 

grid thickness 

.3 cm 

•5 mm 

.3-4 cm 

pairs/Fourier comp. 

1 

2 

1 or 2 

image ambiguity 

no 

no 

yes (w/ 1 pair) 

centerline image cap. 

yes 

yes 

no (w/ 1 pair) 

rotation required 

no 

no 

yes 

time resolution 

1 s 

1 s 

TBD 

gnd separation 

5 m 

1.4 m 

5 m 

spatial frequencies 

2 

6 

12 

spatial resolution 

4 arc sec 

5 arcsec 

2 arcsec 

finest slit width 

.0125 cm 

.0060 cm 

.0050 cm 

detector material 

Nal 

Nal 

Germanium 

detector area/pair 

100 cm 2 

4 cm 2 

TBD 

det. elements/pair 

8 

1 

1 

total sensitive area 

4800 cm 2 

256 cm 2 

TBD 

energy range 

10-100 keV 

15-100 keV 

10 keV -20 MeV 

field of view 

full sun 

full sun 

full sun 


Table 5.1: Our SMC Compared To YOHKOH And HESP 
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Recalling that at large distances photon trajectories from the same point 
are essentially parallel, a photon may fall anywhere randomly on the surface of 
the telescope. Thus, each photon arriving at the face of the telescope has associ- 
ated with it two angles representing the source location in the sky. Each photon 
arrived at a randomly selected grid (48 total) in the top grid plane with a ran- 
domly selected subset of angles (relating it back to the set of angles describing 
the extended source) at a randomly selected local (x,y) coordinate on the grid. 
The model was then tested to see if the photon survived passage through the first 
grid, calculated its (x, y) position on the second one, and tested to see if it sur- 
vived passage through the second grid. This produced a detector response which 
can be reduced to Fourier components and converted to an image as described 

previously (Figure 5.6). 

C. RESULTS 


Figure 5.6 illustrates the solar flare input model and gives a close-up view 
of the synthetic flare configuration chosen for simulation. The flare is made up 
of 21 pixels of equal size. An odd number was chosen so as to provide a point at 
the apex of the loop and one at the footpoints. 

Figure 5.7 illustrates the dirty image of a full flare created from the synthetic 
flare by allowing each of the 21 pixels to have an equal probability of producing 
a photon. Figure 5.8 shows a cleaned version of the flare. Due to the fact that 
only two spatial frequencies were being measured, the telescope seemed to have 
greater difficulty with the extended, more complicated, image than it did with 
the simpler two footpoint test image. As one can see the spatial resolution of the 
instrument is about 4 arc seconds, consistent with the selected spatial frequency. 
Approximately 10 4 photons per 4 by 4 arc second telescope resolution cell are 
required by the telescope to form a meaningful extended image (at the Earth). 
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Figure 5.5: SMC Numerical Model 
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Figure 5.6: SMC And RMC Flare Model 

This model consists of 21 pixels equally spaced to form a half circle. (An odd 
number was chosen so as to provide a pixel at the top of the loop and one at the 
footpoint.) If most or all of the pixels have an equal probability of emitting a 
photon the telescope with 4 arc second resolution will see the loop as an extended 

source. 
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Figure 5.7: SMC Dirty Image Of Loop 

The y-axis is in arc minutes and the x-axis is in hours, minute, seconds from 
Sun center. In this dirty image, there was an equal probability that all parts of 
the loop would emit a photon. An observer seeing this image for the first time 
would note the symmetrical pattern and suspect that cleaning would reveal a 
significant image at the center of the pattern. 
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Figure 5.8: SMC Clean Image Of Loop 

This shows the effects of the cleaning and zooming features in AIPS to better dis- 
play the resultant image. This is the final output of our basic spatial modulation 
collimator and indicates the promising potential of this instrument for viewmg 
complicated shapes. The y-axis is in arc seconds and the x-axis is in hours, min- 
utes and seconds from Sun center. The imaging threshold for extended sources 
such as these was found to be = 10 4 photons per 4 by 4 arc second telescope 

resolution cell (at the Earth). 
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These instruments have an associated point spread function similar to that 
for any optical telescope (or a beam in radio astronomy). Even for the basic 
telescope configuration selected, the point spread function exhibits the sharp 
central peak necessary for high resolution imaging (Figure 5.9 ). 

Next choosing a two footpoint object for use as a standard test configuration, 
the footpoints were displaced from Sun center in ever increasing increments until 
the test image was seriously degraded. In the example shown in Figure 5.10, the 
footpoints show little change at 15 arc minutes as compared to those at 30 arc 
seconds. Beyond 17 arc minutes the image began to seriously degrade; however, 

17 arc minutes (34 arcminute diameter) should be sufficient to provide full Sun 
coverage. 

During the development of the SMC, one systems-related question frequently 
arises when applying the conceptual technology to real structures and pointing 
systems. One of the most pressing issues has been the subject of twist defined as 
the rotation of one grid with respect to the other due to elastic deformation of the 
structure supporting, separating, and serving to keep the grids aligned. Private 
conversations with the designers of HEIDI have indicated that they expect to 
control twist to the sub arc second level. 


Twist was found to become increasingly significant above 2 arc minutes (Fig- 
ure 5.11). For twist less than this, the images were stable. Current applications 
are expected to limit twist to a several arc seconds, 

Fourier telescopes are also limited in dynamic range. In other words, the 
telescope takes the the brightest source in the field of view as a reference baseline 
and images only those other sources within approximately a factor of 10. If a 
source is weaker than a factor of 10 from the brightest source, then it will most 
likely be lost in the artifacts produced by the telescope. Figure 5.12 illustrates 
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Figure 5.9: SMC Point Spread Function (Beam) 

This basically shows how the telescope would image a single point source. Clearly, 
having a bright central peak which stands out above the noise is crucial to suc- 
cessful imaging. In radio astronomy terminology, this would be described as a 

beam. 
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Figure 5.10: SMC Full Sun Field Of View Demonstrated 

This shows two footpoints near the limb of the Sun and verifies that the spatial 
modulation collimator can accomplish full Sun field of view imaging for the ge- 
ometry used. It also illustrates the point that the effective field of view of the 
telescope is much less than the geometrical field of view. The y-axis is in arc 
minutes and the x-axis is in hour, minutes, and seconds from Sim center. 
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Figure 5.11: SMC Twist 

This shows the spatial modulation collimator operating with twist between the 
two grids. The tour dirty images show that little degradation occurs below 2 arc 
minutes. Current designs expect to hold twist to a tew arc seconds. 
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this aspect with a two footpoint minimum flux input model. The brightest 
footpomt’s probability of emitting a photon is 10 times that of the weakest. 

D. SUMMARY AND CONCLUSIONS 

This chapter has described the construction and verification of a simulated 
Fourier telescope for viewing the Sun. No rotation is required for this instrument; 
hence, we call it a spatial modulation collimator. A very basic telescope was 
constructed measuring only two spatial frequencies, yet it did surprisingly well. 
Test highlights axe discussed below. 

The emission profiles from Chapters II and III were also viewed by the 
telescope to determine in general if Fourier telescopes could resolve the differences 
in the two sufficiently. Indeed, the profiles when viewed by the spatial 
modulation collimator did result in significantly different images. This 
will be discussed further in Chapter VII. 

One of the answers that the numerical model provides is the number of 
random photons from a synthesized object required to make an image. The 
telescope simply counts photons (i.e., builds Fourier components statistically) 
to create an image. It leaves to the observer the task of managing bin sizes 
and energies. Once the telescope imaging threshold requirement is known, one 
can compare the threshold to real observations. Generally, for the SMC, = 10 4 
photons per 4 arc second by 4 arc second telescope resolution cell were sufficient 
to adequately image Model T and Model N-T signatures, respectively. Imaging 
above this threshold guaranteed that image location, shape, and orientation were 
reproduced correctly constrained only by the spatial resolution capability of the 
telescope. As the photon level was decreased, shape degraded first, followed 
by orientation, and then finally position. Generally, the threshold marked the 
beginning of a gradual deterioration rather than a sharp break. 
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Only two footpoints are emitting in this test. In the first three images, the 
probability of the upper footpoint emitting a photon was 1, 5, and 10 times 
greater than that of the lower point, respectively. In the fourth image, the 
probability of the upper point emitting a photon is 20 times greater than that 
of the lower point. Clearly, the Fourier telescope is limited in imaging contrasts 

greater than 10:1. 
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The field of view of the telescope was checked. In general, one would prefer to 
have full Sun coverage with the telescope to simplify pointing problems; however, 
this may be categorized as a highly desirable requirement rather than an absolute 
one. In the case of the spatial modulation collimator, full Sun coverage was 
achieved for the 10 x 10 cm grid pairs with 5 meter separations; however, this 

effective field of view was only 25% the geometrical field of view as normally 
calculated for optical instruments. 

The telescope was found to be very tolerant to random noise on the detector. 
This was to be expected since the information necessary for imaging is taken 
from the relative modulated signal rather than the total, absolute one. In other 
words, the information for imaging is found in the modulation which rides on 
top of the somewhat uniform random noise curve. Naturally, systematic noise 
will defeat this instrument as it will most others. The tolerance to noise is also 
necessary to the successful performance of the instrument in that one would 

expect a higher level of random noise to be present due to the large detector area 
required. 

One note of caution is that the telescope generally will try to provide a 
solution even m the absence of sufficient flux. An observer attempting to image 
a flare in low flux conditions may get a low quality image that may be misleading. 
Care should be taken in drawing major conclusions from limited flux images. 
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VI. THE ROTATION MODULATION COLLIMATOR 


A. INTRODUCTION 

In the last chapter we contructed a simulation of a spatial modulation col- 
limator to view the emission profiles from the Model N-T (Chapter II) and the 
Model T (Chapter III). In this chapter, we shall construct a simulation of a ro- 
tating modulation collimator simulation for viewing the same profiles (Appendix 
D). In addition, we shall discuss various telescope limitations and degradations. 

Rotating the telescope allows the design of the SMC discussed in the previous 
chapter to be simplified. First, the grid pairs required now are (n,n) grids to 
modulate the incoming photons as opposed to the SMC (n, n + 1) scheme. 
Second, while we retained the 48 grid pairs (4800 cm 2 ), the slit spacings (0.0125 
and 0.0275 cm), and the grid separation (5 m), the detector for each grid pair 
now only needs to be composed of a single element. By rotating this design at 10 
rpm the (u, v) plane could be mapped every 0.25 seconds. Clearly, this approach 
results in a much simpler engineering design with the trade-off that rotation is 
now required. In other words, the resulting image will be integrated over the 
time required by the telescope to sequentially map the («,u) plane. 

B. THEORY 


Figure 6.1 illustrates a ray trace model for the RMC as well as its basic 
operation. The rays will penetrate all the way through both grids to the detector 
or they will be stopped by one of the grids. Moving the point source with respect 
to the grids will produce a regular sawtooth pattern as shown m Figure 2. 

Now, if one rotates the grids, the detector response may be found by ex- 
amining the vector A which denotes a ray of fight passing through a grid at an 
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Figure 6.1: RMC Ray Trace Model 



This illustrates the idea that a ray coming from any direction (rotating tele- 
scope) may be divided into a contributing component across the slits and a 
noncontributing one along the slits. The strength of the signal on the detector 
at angles (0, <f>) will correspond to the value of the regular sawtooth at the an gle 
(«). 
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Figure 6.2: RMC Sawtooth 


This is a regular sawtooth produced by typical grid geometries. Note the maxima 
and minima are separated by very small angles; hence, only a small movement 
of the source across the sky will produce a sawtooth pattern. 



arbitrary angle. A may be reduced to two components parallel and perpendic- 
ular to the slit direction. The amplitude of the component parallel to the slit 
does not contribute to the detector pattern while the value of the perpendicular 
component is taken from the regular sawtooth curve (Figure 6.2) corresponding 
to an offset angle a from the pointing axis. 

Given the angle of rotation of the telescope about its line of sight <j> and given 
the angle of the source 9 from the pointing axis, the angle o may be calculated. 
From Figure 6.1, we may write the following expression, 

tan a = tan 9 cos <f> ^g 2) 

giving 

a = tan -1 [tan 9 cos $]. (g 2 ) 

For small 9 this reduces to 

a = 9 cos <f>. (g >3 ). 

From Figure 6.3, one observes that for a regular sawtooth the zeroes of the 
intensity are located at 


for 


( a min)n = tan 


1 


+ l)s 


(6.4) 


n = 0,1, 2, 3,.... 

where d is the distant between grid planes and s is the slit width. Similarly, the 
maxima are located at 
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for 


1 / Zdt to i 

(^mai)n = \ d ) 


(6.5) 


n = 0, 1, 2, 3, ... 

Equating these two expressions with equation (6.1) gives the maxima and 
minima of the amplitudes for the corresponding rotational angle of the telescope, 



where 


n = 0, 1, 2, 3, ... 

These expressions allow the computation of the RMC ray trace detector 
response. Figure 6.4 illustrates the compressed sawtooth found from this ray 
trace model. These ray traces, while not used directly for imaging, provide 
useful information regarding the rotating collimator. For the ray trace variable 
sawtooth, as the point source angle from the pointing axis, 6, increases, the 
sawtooth frequency will increase. For a specific 0 it as the angle, <f>, (Figure 6.1) 
is varied from 0 to 2 t r, the location of the frequency minima of the curve in Figure 
6.3 will also vary between 0 and 2 tt in a corresponding manner. In addition, the 

sawtooth repeats itself for the second 180° of rotation. 

At this point we leave the ray trace model and (similar to the SMC) use 
a more realistic random, photon counting model to obtain images. As will be 
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UPPER GRID 


LOWER GRID 


Figure 6.3: Grid Pair Maxima and Minima 

For a given geometry, this illustrates that for a point source moved across the 
sky perpendicular to the slits of a grid pair, sets of associated angles exist cor- 
responding to intensity maxima and minima. Note, changes in the slit width, 
s, and/or the grid separation, d, will cause significant changes in the repetitive 
pattern of maxima and minima seen by the detector. 
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ANGLE OF ROTATION (radians) 


Figure 6.4: RMC Ray Three Compressed Sawtooth 

This variable sawtooth trace obtained from rotating the telescope through one 
full rotation while observing a point source located at an angle t off the line of 
sight axis rotated by a second angle fo in a plane perpendicular to the Kne of 
sight axis. Although somewhat subtle, one can readily « that A, is The 

angle » can only be determined by the level of compression of the sawteeth. As 
the point source is moved further away from the line of sight the compression 
of the sawteeth rapidly exceeds the line width of the printer giving a sohd black 


square. 
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seen later, the ray trace model was valuable in the development of the RMC 
simulation code in verifying that the form of the detector response was correct. 
Figure 6.5 illustrates the numerical model used for the RMC and while it is similar 
to the SMC model, there are distinct differences. For example, the photon is now 
associated with a random rotation angle of the telescope, and two grid pairs are 
used to provide the real and imaginary components for one point in the (u,v) 

plane. Table 6.1 compares this RMC to Solar-A (YOHKOH) and to the proposed 
HESP instrument. 

The detector output for our random, photon counting RMC model is illus- 
trated m Figure 6.6 for both the real (cosine) and imaginary (sine) components 
of a single point source. Each spatial frequency is represented by two curves. 
Each curve represents essentially all Fourier components for that particular spa- 
tial frequency. These curves approximate sinusoids and thus, for a specific time 
of rotation, the real and imaginary values of the respective Fourier components 
can be read directly. In other words, these four curves contain all the information 
necessary to create an unambiguous image of the associated point source. In the 
same manner as the SMC, the RMC via AIPS first provides a dirty image which 
must be processed the obtain a clean , final image (Figure 6.7). 

The telescope point spread function is illustrated in contour in Figure 6.8 
and in a three-dimensional plot in Figure 6.9. The dominant central tendency 
indicates that the telescope will image well. The surrounding higher order max- 
ima are a source of noise in the resultant image. AIPS provides software to 

allow these images to be deconvolved allowing the observer to remove this noise 
numerically. 


C. RESULTS 

Figure 6.10 shows a clean image of a full 21 pixel loop. Note, the image 
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grid material 
grid pairs 
grid thickness 
pairs/Four. comp- 
image ambiguity 
centerline image cap. 
rotation required 
time resolution 
grid separation 
spatial frequencies 
spatial resolution 
finest slit width 
detector material 
detector area/pair 
det. elements/pair 
total sensitive area 
energy range 
field of view 


RMC 

YOHKOH | 

tungsten 

tungsten 1 

48 

64 

.3 cm 

.5 mm 1 

2 

2 

no 

no 

yes 

yes 1 

yes 

no 

1 s(10 rpm) 

1 s 1 

5 m 

1.4 m 1 

2 

6 

4 arcsec 

5 arcsec 1 

.0125 cm 

.0060 cm I 

Nal 

Nal 1 

100 cm 2 

4 cm 2 

1 

1 

4800 cm 2 

256 cm 2 

10-100 keV 

15-100 keV 

full sun 

full sun 


Table 6.1: Our RMC Compared To YOHKOH And HESP 
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Figure 6.5: RMC Simulation Model 
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Figure 6.6: RMC Real And Imaginary Photon Sawtooths 

This shows the real and imaginary traces for a rotating modulation collimator. 
All Fourier components may be obtained from these four curves. These four 
curves in fact constitute the image of a point source located 90 arc seconds from 

the line of sight axis. 


139 



00 00 04 03 02 01 (S3 59 59 SB 57 56 

RIGHT ASCENSION 


Figure 6.7: RMC Cleaned Loop 

An example of a clean image of a loop. In this case, 21 pixels constitute the loop 
each having an equal probability of emitting a photon. Clearly, a basic rotating 
modulation collimator having only two spatial frequencies shows promise for 
imaging complex shapes. The imaging threshold for complex images such as this 
was found to be = 10 4 photons per 4” by 4” telescope resolution cell (at the 
Earth). This performance is similar to that of the SMC discussed in Chapter V. 
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Figure 6.8: RMC Point Spread Function (Beam) 

This image of the rotating modulation collimator’s point spread function was 
included to illustrate the need for other graphical approaches and to point out 
the symmetry provided by the telescope. This symmetry is important to the 


construction of an image. 



Figure 6.9: RMC Beam 3-D Representation 


The bright central peak is necessary for successful imaging in that the point 
spread function illustrates the ability of the telescope to recreate a point source. 


2 


is rotated 90° from Figure 6.7, illustrating the telescope’s ability to image the 
loop rotated at any angle. As for the SMC, a stable RMC image requires Si 
10 4 photons per 4 arc second by 4 arc second telescope resolution cell (at the 
Earth). Again, only about 25% of the available photons are reaching the detector. 
In addition, the telescope successfully accomplished full Sun imaging achieving 
about 4 arc seconds spatial resolution. And, as with the SMC, the RMC is highly 
tolerant of random noise on the detector since the modulation contains all the 

necessary information. 

Figure 6.11 illustrates the twin point source dirty image of the RMC. Note 
the distinct differences between the patterns of the noise around the point sources 
as compared to those in Chapter V. Figure 6.12 illustrates the effects of twist 
on the RMC. Note, image integrity is maintained up to 2 arc minutes of twist. 
Above this limit the image begins to degrade in shape. At 10 arc minutes, the 
image is lost completely. Private conversations with the designers of HEIDI have 
indicated that they expect to control twist to the arc second level. Clearly, given 
this type of alignment control, twist will not be a factor in the imaging quality 

of the telescope. 


D. SUMMARY 

In this chapter, we have constructed and verified an RMC photon counting 
simulation. The design parameters chosen were essentially the same as those of 
the SMC. This RMC employed 48 grid pairs and required only 15° rotation to 
completely map the (u,v) plane (2.5 rpm for 1 second temporal integration). A 
single fourier component required the output from two grid pairs (i.e., a real and 
an imaginary). Low flux imaging performance was found to be approximately 
the same as the SMC. In addition, for high flux levels snapshot images could be 
attained for time intervals on the order of a fraction of a second. 
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RIGHT ASCENSION 


Figure 6.10: RMC Rotated Loop 

In this case, a uniformly emitting loop was rotated to test the telescopes ability 
to image it at an different angle. Clearly, even a basic rotating modulation 
collimator shows promise for imaging complex shapes. 


144 





Figure 6.11s RMC Footpoint Image 

A dirty image of two footpoints is provided here to illustrate the difference m 
noise patterns between the RMC and the SMC. However, both patterns while 
different axe symmetrical in nature consistent with their point spread functions. 
The observer would expect to find a real image at the center of such a pattern. 
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Figure 6.12: RMC Twist 

These four dirty images of two footpoints clearly show little or no significant 
change below 2 arc minutes. Designers believe that twist may be reasonably 
restricted to a few arc seconds in actual designs. 
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It should be noted that another advantage of the RMC approach is to allow 
the use of fewer grids at the expense of reduced sensitive area. For example, four 
grid pairs is the minimum number which could be used to accomplish imaging 
similar to what we have done with the SMC and the RMC (i.e., two spatial fre- 
quencies), all other design parameters remaining the same. A 1 second temporal 
integration requirement would mean that this design would have to be rotated 
at 30 rpm, and the threshold for imaging extended sources would increase by an 

order of magnitude. 

In conclusion, the 48 grid pair RMC was found to be virtually equivalent 
to the SMC in performance. Both telescopes offer advantages and disadvantages 
which can translate to flexibility for the designer. Chapter VII will discuss the 
performance of a Fourier telescope viewing the flare profiles derived in Chapter 

II and III. 
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VII. FOURIER IMAGES OF SYNTHETIC FLARE PROFILES 


The primary objective of this work is to address how well Fourier telescopes 
can image hard x-ray sources on the Sun. In answering this question, two ele- 
mental loop models, the Model T (Chapter II) and the Model N-T (Chapter III), 
both employing typical geometries and x-ray flare characteristics, were investi- 
gated quantitatively. The dependence of the x-ray emission intensity upon time, 
energy, and loop position was determined (Chapters II and III). 

Given the fact that significant spatial differences were found at particular 
energies and times in the two hard x-ray models, the question now becomes how 
well can a Fourier telescope resolve these. A basic SMC (Chapter IV) and a 
basic RMC (Chapter V) employing typical geometries (e.g., grid separations, slit 
widths, detector areas) and individual, random photon counting were constructed 
numerically. AIPS was modified to process the output from these simulations and 
to reconstruct clean images (Chapter III). Both Fourier telescopes with similar 
designs (e.g., sensitive areas, spatial frequencies, etc.) were found to perform 
equally well against test cases. In this chapter, we investigate the performance 
of the Fourier telescope against the Model T and the Model N-T profiles to 
determine if the differences found in those profiles may be imaged. 

Figure 7.1 qualitatively summarizes the hard x-ray emission from the Model 
N-T as a function of time, energy, and loop position. At low photon energies and 
early m time, the predominant emission originates in the footpoints of the loop. 
Later in the event, thermal effects dominate to produce a bright central peak. 
At higher energies, the twin footpoint emission dominates throughout the event. 
At the end of the event at low energies, after beam shutoff, a filled loop remains 
while no significant emission is observed at higher energies. 

Similarly, Figure 7.2 qualitatively summarizes the hard x-ray emission from 
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Photon Energy (keV) 



Figure 7.1: Model N-T Intensity Profiles 

This illustrates qualitatively the intensity as a function of flare position at several 
energies for several different times (Chapter II). The broken lines indicate that 
emission is occurring at that location but that its intensity is below the Fourier 

telescope imaging threshold. 
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the Model T as a function of time, energy, and loop position. At low photon 
energies and early in time, the predominate hard x-ray emission originates in the 
apex of the loop and spreads to adjoining pixels as the hot kernel grows. The 
intensity decreases in time reflecting the rapid cooling of the emitting plasma 
within the kernel. Nonthermal emission is predominant in the footpoints but is 
overshadowed by the thermal emission from the kernel. The triple emission point 
signature is prevalent from 10-40 keV. At higher energies and very early in time, 
a brief pulse of thermal emission is noted at the beginning of the event, but it 
disappears quickly leaving only emission from the footpoints. At higher energies 
and later in time, the Model T mimics the Model N-T to a degree in that weak 
emission from the footpoints is all that is visible. 

Thus, we see that while ambiguities exist between the two models, especially 
at higher energies, there is sufficient information at lower energies (i.e., 10-40 keV) 
to distinguish between the two x-ray models. This information must be in the 
form of sequential observations with approximately 1 second integration times 
and 4 arc second spatial resolution starting preferably just prior to the start of 
and continuing throughout the impulsive phase. 

Figure 7.3 illustrates the Fourier telescope’s view of the Model N-T x-ray 

flare at four different times at 10 keV based upon a 1 second integration time 
and 10 keV bins. 

Figure 7.4 shows the telescope’s view of the Model T at selected times at 10 
keV. Early in time, the telescope will view this model as a single bright source. 
Due to dynamic range limitations, this single source will dominate the image As 
the event proceeds the bright central source spreads and begins to diminish in 
intensity. Since there was only one energy release in this model, the loop rapidly 
cools and imaging is lost after 10 seconds into the event. The four sequential 
images shown in Figure 7.4 axe based upon a 1 second integration time, 10 keV 
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Photon Energy (keV) 



Figure 7.2: Model T Intensity Profiles 

This illustrates qualitatively the intensity as a function of flare position at several 
energies for several different times (Chapter III). The broken lines indicate that 
emission is occurring at that location but that its intensity is below the Fourier 

telescope imaging threshold. 
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ao SEC INTO THE EVENT >>60 SEC INTO THE EVENT 


(AFTER BEAM SHUTOFF) 


Figure 7.3: Model N-T Hard X-ray Images At 10 keV 

These four images were formed from the quantitative profiles for this x-ray model 
and are based on a 1 second integration time, 10 keV bins, and a 4800 cm 2 
sensitive area. 
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bins, and a 4800 cm 2 sensitive area. 

Figure 7.5 shows a comparison of how a Fourier telescope would image the 
Model N-T and the Model T profiles at selected times at 40 keV. At 1 second 
into the event, distinctive images are observed. However, the Model T rapidly 
fades from view. This implies a requirement that imaging take place early in 
the impulsive phase. Again, images are based upon 1 second integration times, 
10 keV bins, and a 4800 cm 2 sensitive area. Clearly, the Fourier telescope 
can potentially observe differences in the signatures from 10 to 40 keV 
early in the event. 

Figure 7.6 shows a comparison of the two models at 100 keV. Only Model 
N-T will be clearly imageable 1 second into the event at this energy, although 
Model T may faintly mimic its signature. 

Clearly, from Figures 7.1 and 7.2 for photon energies at and below 40 keV, 
differences are present in the hard x-ray model profiles. Figures 7.3 and 7.4 
illustrate that these differences also are also present to some degree in the images 
at 10 keV. Early in the event, the Model T is imaged as an expanding bright 
point at the apex of the loop (dynamic range suppresses the footpoints) while 
the Model N-T gives the distinctive twin footpoint signature. 

Briefly at the beginning of the event at 40 keV (Figure 7.4), the Model T 
provides a triple point signature but soon evolves to a state which is nommageable 
while the Model N-T provides a twin footpoint signature throughout the period 
in which a beam is present. Above 50 keV at all times during the event the 
Model T profile generally begins to look like that of the Model N-T in that 
emission from the footpoints is dominant. Any small differences in these profiles 
will be further suppressed in the associated Fourier telescope image. Thus, the 
optimum energies for viewing flares as indicated by the images of hard 
x-ray model profiles will be from 10 to 40 keV. Primarily, this is because 
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Figure 7.5: Model N-T And Model T Hard X-ray Images At 40 keV 

The upper two images are the Model N-T and the lower two are the Model T. 
Images are from the quantitative hard x-ray profiles derived in Chapters II and 
III and are based on 1 second integration time, 10 keV bins, and a 4800 cm 2 


sensitive area. 
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Figure 7.6: Model N-T And Model T Hard X-ray Images At 100 keV 


The Model N-T images show a clear footpoint signature while the Model T images 
seem to weakly mimic the Model N-T. Images are based on 1 second integration 
time, 10 keV bins, and a 4800 cm 2 sensitive area. 
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the thermal component of the hard x-ray emission in each model, m general, is 
crucial to providing discriminating differences in the images. 

Temporal resolution is also seen to be important in comparing the Model 
N-T and Model T profiles. For higher energies 40 keV and above, a few seconds 
should suffice to provide a reasonable image. However, for energies less than 
40 keV especially in the Model T, significant changes in emission signatures 
are occurring on temporal scales of a second or less. Therefore, our choice of 
a one second integration time for our telescopes was a good one. In addition, 
the capability to provide shorter integration times is indicated as a desirable 

requirement. 

In conclusion, Model N-T and Model T each exhibit emission pro- 
files which, given a favorable viewing angle, can generally be imaged 
by a Fourier telescope with reasonable design parameters. 
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VIII. SUMMARY AND CONCLUSIONS 


The primary objective of this work was to address how well Fourier tele- 
scopes can image hard x-ray sources on the Sim. In answering this question, two 
elemental loop models, the Model T (Chapter II) and the Model N-T (Chapter 
III), both employing typical geometries and x-ray flare characteristics, were in- 
vestigated quantitatively. The dependence of the x-ray emission intensity upon 
time, energy, and position was determined (Chapters II and III). The spectra of 
the Model N-T and the Model T were consistent with observations. However, 
no significant differences were found between the two spectra further illustrating 
the need for an imaging approach. 

Given the fact that significant spatial differences were found at particular en- 
ergies and times in the two models (Chapters II, III, and VII), the question then 
became how well a Fourier telescope could resolve these. A basic SMC (Chap- 
ter IV) and a basic RMC (Chapter V) employing typical geometries (e.g., grid 
separations, slit widths, detector areas) and individual, random photon counting 
were constructed numerically. AIPS was modified to process the output from 
these simulations and to reconstruct clean images (Chapter III). Both Fourier 
telescopes with similar designs (e.g., sensitive areas, spatial frequencies, etc.) 
were found to perform equally well. The profile differences discussed previously 
were then provided as input for the Fourier telescope and imaged using AIPS. 
The results of this work were discussed in Chapter VII. In short, the Fourier 
telescope was found to be able to resolve these differences early in the event and 
from 10 to 40 keV. 


A. CONCLUSIONS 

Both basic telescopes measuring only two spatial frequencies required ap- 
proximately 10 4 photons per 4 x 4 arc second telescope resolution cell (at the 
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Earth) to image clearly an extended source. While both types of telescope pro- 
vided different trade-offs, both were viewed as equally promising for viewing the 
Sun at 1 arc second spatially and at 1 second temporally. 

Image degradation due to low flux levels was found to occur in stages. As 
the flux level gradually decreased, image sharpness destabilized first followed by 
image orientation and last by image position. Thresholds were chosen as that 
flux level at which image sharpness first showed significant loss of sharpness or 
resolution in the extended source case. Hence, even for flux levels below the 
threshold some information may be obtainable from the source. 

Rotation is not an absolute requirement in constructing Fourier 
telescopes. The spatial modulation collimator could easily be integrated into a 
three-axis pointing system, thus allowing complementary instrumentation to be 
accommodated. 

Numerical simulations incorporating AIPS are very useful in understanding 
the performance of Fourier telescopes. For example, we discovered that effective 
fields of view for both the SMC and RMC are significantly smaller (i.e., 1:4) 
than the geometrical ones. This means that overall grid areas and hence detec- 
tor areas are strongly coupled for full Sun viewing to boom/cannister lengths. 
The individual detector sizes for both the SMC and the RMC were 10 cm and 
provided a 34 arc minute field of view with 4 arc seconds of spatial resolution. 
The grid separation used was 5 meters. 

The Fourier telescope offers the observer in principle a large bandpass. In- 
deed, HESP is planned to operate from 10 keV out to 4 MeV. In practice, the 
Fourier telescope’s applicability is limited in both the long and short wavelength 
regimes. In the former case, it is diffraction limited and in the latter by detector 

and grid limitations. 

The trade-off to the wide bandpass using these types of telescope is that 
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the dynamic range is inherently limited (Chapter V). In other words, in the 
relative manner in which the telescope ignores random noise, it also sees only the 
brightest points within its field of view to within about an order of magnitude in 
brightness. For an extremely bright central peak with bright footpoints against 
a background almost as bright, the telescope will image the bright central peak 
and the footpomts will be lost in the artificial noise generated due to the limited 
coverage of the (ti, v) plane. There are some techniques which may be developed 
in future work which may improve upon this basic limitation. 

Detailed numerical modeling is a necessary step toward success- 
fully optimizing the design of these instruments. No purely analytical 
approach is gomg to be able to determine the number of photons which pro- 
vide a good image simply because a good image is determined by the observer 
comparing the image to the original synthetic object. Also, through modeling 
the (u, v ) plane coverage can be tested and an optimum verified. Another good 
example is twist. Only an observer looking at images illustrating the effects of 
different levels of twist can decide at which point twist becomes unacceptable, 
and only a numerical simulation allows one to look at twist independent of the 
other degrading effects. 

Hence, both numerical and analytical approachs are useful in understanding 

these instruments. However, the numerical approach provides the proof of the 
pudding. 

Finally , additional considerations for high time and spatial resolution is the 
fact that there may be other hard x-ray models yet to be developed which may 
also be consistent with observations and many loops may be involved in a single 
flare. These loops may be triggered in a somewhat sequential mann er; thus, the 
telescope must be able to have sufficient resolution to resolve these overlapping 
sequences. The proposed HESP mission offers the technology to do just that. 
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B. OPTIMIZATION RECOMMENDATIONS FOR THE FUTURE 

Different approaches may be taken to optimizing a design of a future mis- 
sion depending upon the cost and schedule contracts involved. For an orbital 
platform, one would prefer to be able to look not only at hard x-rays but also at 
many other wavelengths and spectral lines and perhaps have a colocated mag- 
netograph all operating simultaneously. Generally, one would want three-axis 
pointing since that would simplify the design of most of the instruments. This is 
a high cost approach similar to the Hubble Space Telescope requiring a decade 
or more of development time. 

For this case, a spatial modulation collimator would be the optimum choice 
for the hard x-ray imager. Within reasonable geometries, six spatial frequencies 
representing 1, 3, 6, 13, 27, and 56 arc second resolution would be measured 
using 72 grid pairs each measuring 10 x 10 cm 2 for a total collecting area of 
7200 cm 2 (less than a 1 meter diameter aperture). Each spatial frequency would 
be measured by 12 grid pairs oriented from 0° to 180° in 15° steps. 

The grid separation would be 5 meters. The grid material would be tungsten 
and an aluminum filter would shield the aperture to prevent detector saturation 
from low energy flux. High purity germanium would be the best choice for the 
detector which would be composed of an array of seven elements. Telemetry 
requirements would be significantly reduced as each transmitted image would 
simply be an array of 2 X 72 numbers, the real and imaginary components of the 
brightness distribution. The hard x-ray telescope would image the Sun from 10 
keV to as high an energy as cost and schedule constraints would allow. 

For the moderate cost case, a dedicated Fourier Telescope is recommended. 
This would be an RMC to minimize the cost involved in providing the detector 
arrays as required by the SMC. In this case, 48 grid pairs would provide imaging 
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at spatial frequencies represented by 1, 3, 6, 13, 27, and 56 arc seconds. Each 
spatial frequency would be measured by 10 grid pairs, 5 real and 5 imaginary, 
oriented from 0° to 180° in 45° steps. The grid sizes would again be 100 cm 2 
separated by a spacing of 5 meters. However, in this approach, only single 
detectors would be required but telemetry requirements would be increased as 
each time tagged photon event would need to be transmitted. The hard x-ray 
telescope would image the Sun from 10 to 100 keV. 

The large detector area per grid pair would ensure full Sun coverage while the 
5 meter separation between grid planes would allow us to use grids with broader 
slit widths, thus reducing the manufacturing risk somewhat. The telescope would 

be rotated at 15 rpm to ensure a 1 second integration time with a snapshot every 
0,5 seconds. 


C. RECOMMENDATIONS FOR FUTURE RESEARCH 


The limited dynamic range of these telescopes may possibly be improved by 
masking the brightest sources and reconstructing the map. In effect, one could 
strip away one layer at a time in an attempt to reveal more diffuse emission 
regions. This would be an important investigation with implications for Solar-A 
and HESP data analysis. 

HEIDI and especially Solar-A should be useful experiments for determining 
the practical applicability of Fourier telescopes to viewing solar flares. In ad- 
dition, other approaches to the hard x-ray imaging problem are being pursued. 
These include the Fresnel zone plate, the Fresnel spiral, the grazing incident mul- 
tilayer telescope (Walker et al., 1988), and a multi-crystal Bragg diffraction lens 
(Smithers, 1991). 
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The Fresnel zone plate and the Fresnel spiral axe also variations of a Fourier 
telescope and investigations have begun and will continue in the future to deter- 
mine the effectiveness of these instruments. An example of the Fourier spiral is 
shown in Figures 8.1 and 8.2. Figure 8.1 is transparent and by superimposing 
it on Figure 8.2, a Moire fringe pattern is produced. Offsetting one spiral with 
respect to the other is equivalent to a single point source illuminating a spiral 
pair at the same angle from the telescope pointing center. 

Initial results indicate that both instruments provide an annular field of 
view. If the designer attempts to enlarge the outer boundary of the field of the 
view, the inner boundary also expands and vice versa. As the point source moves 
further away from the pointing center of the telescope, the fringes become very 
fine and at a small angle from the center the fringes disappear altogether, thus 
d efinin g the annular field of view. Also, the outer boundary falls far short of full 
sun coverage for acceptable inner boundary geometries. 

Figure 8.3 illustrates the Fourier transform of four point sources and Figure 
8.4 illustrates the spiral pair detector pattern when illuminated by four point 
sources. The two dimensional Fourier transform of four point sources produces 
a checkerboard pattern in the (u, v) plane. Clearly, the basic Moire pattern in 
Figure 8.4 suggests the correspondence of the Moire fringes to a Fourier trans- 

form. 

The basic spiral pair response is a Fourier transform superimposed upon 
a Fresnel kernel. Future work should continue to look at this device to see if 
actual images can be reconstructed. The Fourier spiral telescope would probably 
require a minimum of four spiral pairs each requiring a two dimensional array 

detector. 

Hard x-ray detector research for application in Fourier telescope should be 
initiated. The Marshall Space Flight Center (MSFC) offers exceptional facilities 
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Figure 8.2: The Fourier Spiral (Opaque) 

Follow the instructions given for the preceding figure to form the detector pattern 
for a single point source. 
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Figure 8.3: The Fourier Transform Of Four Point Sources 


The two dimensional Fourier transform of four point sources produces a checker- 
board pattern (Bracewell, 1965). 
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Figure 8.4: Four Fresnel Spiral Equivalent Point Sources 

This basic Moire pattern corresponds to a Fourier transform of the four point 
sources shown in Figure 8.3. 
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for this type of research including the AXAF test facility. Innovative means of 
manufacturing and aligning Fourier grids should be investigated. Integrating 
numerical modeling with the flexible boom technology found at MSFC would 
be an interesting and significant experiment. Developing a numerical model of 
Solar-A would allow insight into the possible pitfalls of analyzing the real data. 

Finally, the lunar P/OF Observatory concept with NASA offers the chance 
to place an observatory on the Moon, thus taking advantage of being ground- 
based without an intervening atmosphere degrading or, in the case of x-rays, com- 
pletely absorbing the radiation. This observatory would provide sub arc second 
spatial resolution and tremendously improved sensitivity for Fourier telescopes 
and other instruments. This would allow, for example, weak, small flares to be 
fully resolved both spatially and temporally. Designing the telescope for such 
an observatory should begin with numerical modeling. Other parallel planning 
should be done in the engineering areas such as pointing and control, structures, 
and foundations in the lunar ground. 

Hard x-ray imaging of the Sim is needed and is doable now and should 
provide a plethora of new questions for astrophysicists for generations to come. 
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APPENDIX A 


Model N-T Code 

This appendix contains an example of the code used in Chapter II for calcu- 
lating the hard x-ray emission as a function of time, energy, and loop position and 
is based on one developed by Peng Li at the University of Alabama in Huntsville. 


169 


o * 


C HARD X-RAY SPATIAL STRUCTURE 
C This program is used to calculate the hard X-RAY spatial structure 
c of the non- thermal electron heated solar atmosphere .This part is the 
main program 

C234567891 12345678921234567893123456789412345678951234567896123456789712 

************************************************************************* 
program main 


C Declaration 

INTEGER IJX.SMO 

REALXA(101),Y(101),Zcc(lE3)X(101,lE3).L(lE3).Nm(lE3),T(lE3) 

REAL Jn( 101, 1E3) Jn( 101, 1E3 ) Ep( 101 ) XSlEJEcELUXJnmin 
REALUn,Uth.UoJIN(500)JT(500)^l(lE3)^(lE3)te,e,p 
REAL ITH( 101,1E3)£X( 101.1E3 )Jthto( 101 )4nto( 101 ),TIMENO,g 
REAL JND( 101, lE3),Ithb(500)Jnb(500)Jtha(500)4na(500) 

REAL SLOPJ.1 (500)4N1 ( 101,25),ITH1 ( 101,25),TOT( 101,25) 
REALPOS(20)XLEN(20) 

DATAR/8.52E8/,C/1.38E-23/,e/1.6E-16/,p/l.E38/ 

DATAB/4./XS1/.05/ 

C FIRST DEFINE THE PARAMETERS FOR CALCULATION OF NON-THERMAL CAS 

C TYPE*, ’ENTER THE UN,M,Uo,Flui,Ec’ 

C Accept*, Un.M.Uo.Flux.Ec 
C 

TYPE *, ’REMEMBER TO DELETE ALL NT*DAT;* FIRST’ 

TYPE *, ’ENTER THE UNIT NUMBER:11, 12, 13,14,OR 15’ 

ACCEPT *,UNUM 

OPEN ( UNIT =11 EILE= T 1 ODAT’.STATUS - ’OLD’ ) 

OPEN ( UNIT= 1 2.FILE- T20DAT’, STATUS- ’OLD’ ) 

OPEN ( UNIT- 1 3 .FILE- T30DAT’,STATUS-’OLD’ ) 

OPEN ( UNIT - 14FILE-T40DAT’,STATUS-’OLD’ ) 

OPEN ( UNIT - 15,FILE-’T60DAT’,STATUS- ’OLD’ ) 

C 

OPEN (UNIT- 16EILE=’MODNTl DAT’.STATUS-’NEW’ ) 

OPEN ( UNIT - 1 7JTLE- 'MODNT2 DAT’ .STATU S - ’NEW’ ) 

OPEN ( UNIT= 1 8 FILE- ’MODNT3DAT’,STATUS = ’NEW’ ) 

OPEN (UNIT- 19JFILE- MODNT4 DAT’.STATUS-’NEW’ ) 

C 

OPEN ( UNIT -20.FILE- ’NTl DAT’.STATUS-’NEW’ ) 
OPEN(UNIT=21J 7 ILE=’NT2DAT’,STATUS-’NEW’ ) 
OPEN(UNIT=22FILE=’NT3DAT’,STATUS-’NEW’) 
OPEN(UNIT-23EILE=’NT4DAT’,STATUS- r NEW’) 

OPEN (UNIT-24J 7 ILE-’NT5DAT , ,STATUS-’NEW’) 

OPEN ( UNIT - 25 .FILE - 'NT 6DAT’ .STATU S - ’NEW’ ) 


C 
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o n o o n ~ o o ^ nnoo 


M=1 
EC- 15.0 
FLUX-1.0E11 
UO-11 


DO5J-1.100 
Ep(J)-J 
Y(J)-Ep(J)/Ec 
IF(Y(I).GE.l) THEN 
Y(J)=1 
END IF 
CONTINUE 

READ IN THE DATA FOR NONTHERMAL CASE 
DO 101=1,450 

READ(UNUM,*)Zcc(I),Nm(I),L(I),T(I) 

TYPE*,ZCC(I).NM(I).L(I),T(I) 

continue 


NOW, GET THE COLUMMN DENSITY 

DO 30 1-1,449 
N1(1)=0.000000 

Nl(I+l)=Nl(I)+Nm(450-I)*L(450-I) 

30 continue 

c NOW, CHANGE THE ARRAY ORDER 

DO 35 1=1,449 
N(I)=N1(451-I) 

C IF (LLT.50)THEN 

WRITE( 18 *)TOS=’,I,’COLUMN DENSITY -’J4(I) 
C ENDIF 
35 CONTINUE 


k=2.6e-18 

C 

DO 3 8 J- 1,100 
no 401-1 449 

X(J4)-(2.*K*N(I))/Ep(J)**2 
C IF(JJEQ.10ANDJiE50)THEN 
C WRITE(17,*) , X= , ^(J4) 

C ENDIF 
40 CONTINUE 
38 CONTINUE 
C TYPE*,X(J,I) 

C WE KNOW X(JJ) AND "Ep(J)/Ec"3, WE ARE READ TO USE THE SUBROUTINE 
C 
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C23456789112345678921234567893123456789412345678951234567896123456789712 


CALL INS(X,Y^CS13J1NJN) 

C 

c 

IF (UNUM3Q.1 1 ) THEN 
TIMENO-IO. 

ENDS 7 

IF (UNUM.EQ.12) THEN 
TIMENO-20. 

ENDS 7 

B 7 (UNUM3Q.13) THEN 
TIMENO-30. 

ENDS 7 

IF (UNUM.EQ.14) THEN 
TIMENO-40. 

ENDS 7 

IF (UNUM.EQ.15) THEN 
TIMENO-60. 

ENDS 7 

C 

IF(TIMENOLE.30) THEN 
G-TIMENO/30 
ELSE 

O-2-TIMENO/30 

ENDS 7 

C**************************************************** 

C FOR UNIFORM RESPONSE CASE ONLY 
C 0-1.0 

c***************«****«************«************* 

c 

DO J= 1,100 

c 

c 

DO 50 1-1,449 

JT(J)-285*4/(B+2)/(B-l)*Ftax*G*(Ec**(B-2))*(Ep(J)**(l-B)) 
JND(JJ)-ABS(JN(JJ)-JN(J,I+1) ) 

In(J4)-JND(JJ)*JT(J) 

C TYPE *JN(JJ) 

C IF(J£Q.10ANDJLE50)THEN 
C WRITE( 18,*)’KEV-’J,POS-’,I,’JDIFF-’JND(J,I) 

C ENDS 7 
50 continue 
C 
C 

C NOW IS THE TIME TO CALCULATE THE THERMAL CASE 

C 

C 

DO 70 1-1,449 

EX(J4)-(Ep(J)*e)/(C*T(I» 

S 7 (EX(JJ).GT.87) THEN 

Ith(JJ) — ( 1.2e-l 1 )*Nm(I) **2*L(I) /p/ (Ep(J)*SQRT(T(I) ) ) 
else 

Ith(JJ)-(1.2e-ll)*Nm(I)**2*L(I)/EXP(Ep(J)*e/(C*T(I)))/ 
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1 (Ep(J)*SQRT(T(I) )) 

C 

END IF 
C 

70 continue 

NOW CALCULATE THE EMISSION AS A FUNCTION OF FOOTPOINT PIXEL POSITI 
FOR ELEVEN PIXELS 


DO 1=1,450 
L1(I)=ZCC(I)/1.4E9 

C TYPE *, , I=M,L(I)= , 1.(I),’L1(I)= , JL1(I) 
IF(L1(I)LE.(.096) )THEN 
IN1(J,1)=IN1(J,1)+IN(J4) 
ITH1(J,1)=ITH1(J,1)+ITH(J4) 
TOT(J,l)=INl(J,l)+rTHl(J,l) 

POS(l) =1*1.0 
XLEN(1)=ZCC(I) 

C TYPE*,’I=\I. , ZCC(I)=’,ZCC(I) 

ENDIF 

IF(L1 (I).OT.( .096) -ANDL1 (I) LE.( .191 ) )THEN 
IN1(J,2)=IN1(J,2)+IN(J,I) 

C TYPE *,IN1(J,2) 

rrai(j,2)=iTHi(U)+rra(j,i) 

TOT(J,2)=INl(J,2)+ITHl(J,2) 

XLEN(2)=ZCC(I) 

POS(2)=I*1.0 

ENDIF 

IF(L1(I).GT.(.191)AND.L1(I)LE.(.286))THEN 

IN1(J,3)=IN1(J,3)+IN(JJ) 

rrHi(j,3)=rmi(j,3)+iTH(j,i) 

TOT(J,3)=INl(J,3)+ITHl(J,3) 

POS(3)=I*1.0 

XLEN(3)=ZCC(I) 

ENDIF 

IF(L1 (I) .GT.( .286) ANDi.1 (I) LE.( .381 ) )THEN 
IN1(J,4)=IN1(J,4)+IN(J4) 
ITH1(J,4)=ITH1(J,4)+ITH(J4) 
TOT(J,4)=INl(J,4)+ITHl(J,4) 

POS(4)=I*1.0 

XLEN(4)=ZCC(I) 

ENDIF 

IF(L1 (I) .GT.( .381 ) .AND.L1 (I)LE.( .476) )THEN 
IN1(J^)=IN1(J^)+IN(J4) 

ITH1 ( J5 ) =ITH1 ( J.5 ) +ITH ( J.I ) 
TOT(J,5)-INl(J,5)+ITHl(J,5) 

POS(5)=I*1.0 

XLEN(5)=ZCC(I) 

ENDIF 

IF(L1 (I) .GT.( .476) -ANDil (I)LE.( 571 ) )THEN 
INI ( J,6) =IN1 ( J,6) +IN ( JJ) 
ITH1(J,6)=ITH1(J,6)+ITH(J,I) 
TOT(J,6)=INl(I,6)+ITHl(J,6) 
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POS(6)-I*1.0 

XLEN(6)-ZCC(I) 

ENDIF 

EF(L1(I).0T.(571 ) ANDil(I)UE.(.666) )THEN 
IN1(J,7)-IN1(J,7)+IN(J4) 

ITH1 (J,7)=ITH1 ( J,7)+ITH(J4) 

TOT( J.7) -INI ( J,7) +ITH1 ( J,7) 

POS(7)-I*1.0 

XLEN(7)=ZCC(I) 

ENDIF 

IF (LI (I ) .GT. ( .666 ) ANDil (I ) UE. ( .761 ) )THEN 
IN1(J,8)=IN1(J,8)+IN(J4) 

ITH1 ( J,8 ) -ITH1 ( J,8 ) +ITH( JJ) 

TOT (J,8)-IN1(J,8) +ITH 1 ( J,8 ) 

POS(8)«I*1.0 

XLEN(8)=ZCC(I) 

ENDIF 

IF(Ll(I).GT.(.761)-ANDil(I)lE.(-856))THEN 

IN1(J,9)=IN1(J,9)+IN(J4) 

ITHl(J,9)-ITHl(J,9)+nH(J4) 

TOT(J,9)-INl(J,9)+ITHl(J,9) 

POS(9)-I*1.0 

XLEN(9)-ZCC(I) 

ENDIF 

IF(Ll(I).GT.(.856)-AND4.1(I)iE.(-951))THEN 
INI (J, 10) -INI ( J,10) +IN( J4) 

ITH1 (J.10) -ITH1 (J,10)+ITH(J4) 
TOT(J,10)-IN1(J,10)+ITH1(J,10) 
POS(10)-I*1.0 
XLEN( lO)-ZCC(I) 

ENDIF 

IF(Ll(I).GT.(.951)JVND4-l(I)iE.(1.0))THEN 

IN1(J,11)-IN1(J ( 11)+2*IN(J4) 

ITH1 (J,l 1 )-ITHl ( J,1 1 )+2*ITH(J4) 

TOT(J,l 1 ) -INI ( J,1 1 ) +ITH1 ( J.l 1 ) 

POS(11)-I*1.0 

XLEN(11)=ZCC(I) 

ENDIF 

ENDDO 

C 


IN1(J,12)-IN1(J,10) 

ITH1(J,12)=ITH1(J,10) 

TOT(J,12)=TOT(J,10) 

IN1(J,13)-IN1(J,9) 

ITH1(J.13)-ITH1(J,9) 

TOT(J,13)-TOT(J,9) 

INI (J, 14) -INI ( J,8 ) 

rmi(j,i4)-rmi(j,8) 

TOT(J,14)-TOT(J,8) 

IN1(J,15)-IN1(J,7) 

ITHl(J,15)-rTHl(J,7) 

TOT(J,15)-TOT(J,7) 

IN1(J,16)-IN1(J,6) 
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rmi(j,i6)=rmi(j,6) 

TOT(J,16)=TOT(J,6) 

IN1(J.17)=IN1(J,5) 

ITH1(J,17)-ITH1(J,5) 

TOT(J,17)=TOT(J,5) 

IN1(J,18)=IN1(J,4) 

ITH1(J,18)-ITH1(J,4) 

TOT(J.18)-TOT(J,4) 

IN1(J,19)-IN1(J,3) 

ITH1(J,19)=ITH1(J,3) 

TOT(J,19)-TOT(J,3) 

IN1(J,20)«IN1(J,2) 

ITH1(J,20)=ITH1(J,2) 

TOT(J,20)=TOT(J,2) 

IN1(J,21)=IN1(J,1) 

ITH1 (J,21 )=HH1 (J,l ) 

TOT(I,21)=TOT(J,l) 


DO I- 1,21 

TYPE*^A(I),ITH(JJ),IN(J4),’KEV=’J, , CM=’J 

RE=1.496E13 

ITH1(J4)=ITH1 ( J4) * (4.444E-11 ) 

INI ( JJ) -INI ( J4) * (4.444E- 1 1 ) 

TOT( J4) =TOT( J,I ) * (4.444E-1 1 ) 

BF(J£Q.10)THEN 

WRITE(16*)’KEV=’J,” 

WRITE ( 16* ) TOSITION-M 
WRITE! 16,*)’XA=\XA(I) 

WRITE( 16* )THER PHOT-MTH(JJ) 

WRITE! 16* ) ’NONTHERPHOT-’4N( J,I) 

ENDIF 

ENDDO 

ENDDO 

DOKA-1,11 

WRITE! 17,*)KA,’ ’,’POS=’POS(KA),’ \’XLEN-’,XLEN(KA) 
ENDDO 


DOM-1,10 

IA-M*10 

XINTOT-O.O 

XITHTOT-O.O 

DOKl-1,21 

WRITE! 16* )IA4CUN1 (IAX1 ) JTH1 (IAJCl ),TOT(IA4Cl) 

XINTOT-XINTOT+IN1 (LAJC1 ) / (4.444E-11 ) 
XITHTOT-XITHTOT+ITHl (IAJCl ) / (4.444E-1 1 ) 
ENDDO 

XTOT-XINTOT +XITHTOT 
C WRITE! 17,*)’KEV=’,IA,’XINTOT-’40NTOT 
write! 19* ) , kev=ma;xintot=’4cintot, 

1 ’XirHTOT-’^QTHTOT.’XTOT - ’.XTOT 

ENDDO 
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c 

DOK1-U1 
WRITE(20,*)K1 
DOM- 1,10 
IA-M*10 

IF(IA£Q.10)THEN 
C WRITE(21,*)TOT(IA,Kl) 

C WRITE(21,*)ITH1(IAJC1) 

WRTIE(21,*)IN1(IAJU) 
ENDIF 

IF(IA£Q.40)THEN 
C WRITE ( 22,*) TOT(IA,Kl) 

C WRITE(22,*)ITH1(IAJC1) 

WRTIE(22,*)IN1(IA,K1) 
ENDIF 

IF(IA_EQ.70)THEN 
C WRITE(23,*)TOT(IAJCl) 

C WRITE(23,*)ITH1(IAJK1) 

WRITE(23,* )IN1 (IAJQ ) 
ENDIF 

IF(IA£Q.100)THEN 
C WRTIE(24*)TOT(IAJCl) 

C WRITE(24,*)ITH1(IAJK1) 

WRITE(24 ,* )IN1 (IA.K1 ) 
ENDIF 
ENDDO 
ENDDO 


C234567891 12345678921234567893123456789412345678951234567896123456789712 
C 

DO 45 J-1,100 

Ithb(J)-0 

Inb(J)-0 

Itha(J)-0 

Ina(J)-0 

Ithto(J)-0 

Into(J)-0 

C NOW TO CACULATE THE FOOTPOINT EMISSION OF ONE ARC SECOND 

DO80I-1.M 

Ithb(J)-Ithb(J)+Ith(JJ) 

Inb(J)«Inb(J)+In(JJ) 

80 CONTINUE 

C NOW CACULATE THE TOTAL EMSSION OF T & N FROM THE WHOLE LOOP 

Into(J)-JT(J) 

DO 90 1-1,449 
Ithto(J)-Ithto(J)+Ith(JJ) 

90 CONTINUE 
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C NOW THE TOP OF THE LOOP OF ONE ARC SECOND 


DO 1001=445,449 
Itha(J)=Itha(J)+Ith(J4) 

C IF(JJLE.SMO) THEN 
Ina(J)=Ina(J)+In(J,I) 

C ELSE 

C SLOP=ALOG(ALOG10(Ina(SMO) ) / ALOGlO(Inmin) ) / ALOG(50./SMO) 

C Ina(J)=10**(alogl0(Inmin)*(50./J)**SLOP) 

C END IF 

100 CONTINUE 

C DF(Ina(J).GE.Ina(J-l)) THEN 

C Ina(J)=Ina(SMO)-(Ina(SMO)-Inmin)/(50-SMO)*(J-SMO) 

C 

C234567891 12345678921234567893123456789412345678951234567896123456789712 

C 

C WRITE(Uo,102)Ep(J)4thb(J)4nb(J)4tha(J)4na(J)Jthto(J) 

C WRITE(6*)Ep(J)4thto(J)4nto(J) 

102 FORMAT* 1X.F8.2, IX, 6(E10.5,1X)) 

45 CONTINUE 
CLOSE(UNIT=ll) 

CLOSE (UNIT= 12) 

CLOSE( UNIT-13) 

CLOSE (UNIT- 14) 

CLOSE (UNIT- 15) 

CLOSE* UNIT= 16) 

CLOSE (UNIT- 17) 

CLOSE (UNIT-18) 

CLOSE (UNIT- 19) 

CLOSE (UNIT- 20) 

CLOSE (UNIT-21) 

CLOSE (UNIT- 22) 

CLOSE (UNIT-23) 

CLOSE (UNIT- 24) 

CLOSE (UNIT-25) 


STOP 

END 


C 

SUBROUTINE INS(X,Y,XS1,B,JIN,JN) 
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************************************************************************ 
C HARDX2POR; IS THE PROGRAM CONSIDER THE LOW PART OF X ( J,I ) 

C AND IN ORDER TO GET AWAY THE WIGGLES 
C 

************************************************************************ 

C 

C234567891 12345678921234567893123456789412345678951234567896123456789712 


C THIS PROGRAM IS USED TO CALCULATE THE HARD X-RAY SPATIAL STRUCT 
C WITH THE ENERGY SPECTRUM OF A POWER LAW WITHOUT CUT-OFF. 

C DECLARE 
INTEGER IJN 

REAL JIN( 101 ),Y( 101 )A-( 101 )3,XS( 101) JL( 101,1E3),XS1 
REAL J5( 101, 1E3)J0( 101.1E3) JN( 101.1E3) J1 ( 101,1E3),I2( 101, 1E3) 

REAL Jl( 101, 1E3)J2( 101, 1E3)J3( 101, 1E3)J4(101,1E3) 

C DEFINE THE PARAMETER 
C TP YE*, ENTER YJB,XMAX,XMIN,N,XS1’ 

C ACCEPT*3,XMAX,XMIN,N,XS1 


C 

DO5J-1.100 

A(J)-Y(J)**2 

XS(J)-1/A(J)-1 

JIN(J)-=(4*Y(J)-A(J)**2-3)/12+l/(2-B)-Y(J)/(l-B) 

5 CONTINUE 

C WE NOW READ TO CACULATE THE CUMULATED INTENSITY JN,WE WILL 
C HAVE ONE LOOP AND TWO BLOCKS WITH IF 

DO 8 J- 1,100 
DO 101-1,449 
C TYPE *,’J-’J,I-’J 

C X(U)=10**(XMIN+(XMAX(J)-XMIN)*<I-1)/N) 

C TYPE *,’X-’K(JJ) 

C NOW WE BEGIN THE FIRST IF BLOCK 
C 

IF(X(U)EE_XS(J) )THEN 
C 

C23456789112345678921234567893123456789412345678951234567896123456789712 

C 

J1(U)-A(J)**2/12*(-3*(1+X(U))**2+4*(X(U)+1)**0-1) 

J2(U)-^5*(((1-A(J)*X(U))**15-A(J)**2*(1+X(U))**J>) 

1 +5*A(J)*X(J4) 

1 *(SQRT(1-A(J)*X(U))-A(J)*SQRT(1+X(J4))) 

1 +J*(A(J)*X(JJ))**2 

1 *ALOG(Y(J)*(l+SQRT( 1+X(JJ) ) ) / ( 1+SQRT( 1-A(J)*X(JJ) ) ) ) 

1 +A(J)**2*(1+X(JJ))**2-1) 


IF(X(JJ)JLTJCS1)THEN 
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c 


J0(JJ)=J1(J4)+I2(J4) 
JK(JJ)=J0(J4) /JDST(J) 


C 

ELSE IF(X(J,I) .GE.XS1 ANDX( J4) LE XS ( J) )THEN 
C 

c 

c 

IF(BEQ.4)THEN 

Il(J4)=2*(l-(l-A(J)*X(J4))**L5)/3 
ELSE IF (B£Q.5 )THEN 

I1(J4)«.25*(ASIN((A(J)*X(J4))**.5)-SQRT(A(J)*X(J4) 

1 *( l-A(J) *X( J4) ) ) * ( 1-2*A(J)*X(J4) ) ) 

ELSE IF (BJEQ.6) THEN 

Il(J4)-2*(2-(3*A(J)*X(J4)+2)*(l-A(J)*X(J4))**1.5)/15 
END IF 
C 
C 

J3(JJ)-.5*((A(J)*X(J,I))**(l-5*B)*Il(J,I)-2/(B-2)) 

J0(JJ)=JI(J,I)+J2(J,I)+J3(J,I) 

JN(J,I)-JO(J,I)/JIN(J) 

END IF 

C NOW WE WILL BEGIN THE SECOND BLOCK WITH THE FIRST IF 
C 

c 

ELSE 

C 

c 

IF(B£Q.4)THEN 

I2(JJ)=2*(l-l/(l+X(J4))**15)/3 
ELSE IF (B-EQ.5) THEN 

I2(J,I)-.25*( ASIN(SQRT(X(JJ) / ( 1+X(J4) ) ) ) 

1 -SQRT(X(J4) )*( 1-X(I,I) ) / ( 1+X(J,I) )**2) 

ELSE IF(B£Q.6) THEN 

I2(J,I)-2*(2-(5*X(JJ)+2)/(l+X(JJ))**23)/15 
END IF 
C 
C 

C23456789112345678921234567893123456789412345678951234567896123456789712 

C 

J4(J4)-JIN(J)+Y(J)/(1-B)*(A(J)*(1+X(JJ)))**(3*(1-B)) 

1 -(A(J)*(1+X(J4) ) )**( 1-5*B) / (2-B) 

J5(JJ)-5*(A(I)*X(JJ))**(1-5*B)*I2(J4)-(A(J)*(1+X(J4))) 

1 **( 1-.5*B) / (B-2) 

J0(W)=J4(JJ)+J5(J4) 

IN(J4)=JO(J4)/JIN(J) 

C 

C TYPE VJ4=’44(J4),’J5=’45(J4) 

C TYPE*,’J0=*40(J4),’JN= , 4N(J4) 

C 

c 

END IF 
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c 

IF(J£Q.10ANDJXEJ0)THEN 
C WRnE(18*)’X=’.X(JJ), , JN-’4N(J4),’JO-’JO(J4), 
C 1 , JIN=’4IN(J),’J=’J,’I=’4, , J4=’J4(J ( I), 

C 1 ’J5-’J5(J4) 

ENDIF 

C 101 FORMAT(2X4 J 12.6,2X4 J 12.6) 

10 CONTINUE 
CONTINUE 


RETURN 

END 
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APPENDIX B 


Model T Code 

This appendix contains an example of the codes used in Chapter III for cal- 
culating the hard x-ray emission as a function of time, energy, and loop position. 
In addition, it contains the code used to calculate the evolution of L(t)^ T(t)^ V \ , 
and vico as functions of time. 
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PROGRAM MODT2 


C 

C*********************************************************** 

C THIS PROGRAM CALCULATES THE KERNAL THERMAL EMISSION AND 
C THE FOOTPOINT NONTHERMAL EMISSION FOR THE MODEL T FLARE 
C*********************************************************** 

c 

c*********************************************************** 

C INITIALIZATION 

£********************************************«************** 

C 

REALM XI( 100)ZK 100),SIX( 100) 

REALM WI( 100) ,RI( 100) ,WIT( 100) 

REALM TI( 100), YI(100),XXPOS(20) 

REALM XIT( 100),SI( 100) JCXLEN(50) 

REALM VI( 100, 100),UI( 100,100) 

REALM BI( 100,100), AIT( 100,100) 

REALM XNLEN(25),XNCT(25),XXTE(25) 
REALMXNCP(25.15),XJ(25,15),XXEC(15) 

REALM XXN ( 25, 15 ) .XXI ( 25,25 ) ,XJEL( 1001,501 ) 

REALM XXEP(15),XXIN(25,25),XXJ(25,15) 

REALM XJMX( 15) 

C 

OPEN (UNIT-1 1.FILE-MODT21DAT’, STATUS- NEW) 

OPEN ( UNIT- 1 2 .FILE- MODT22DAT\STATUS- NEW ) 

OPEN ( UNIT- 1 3,FILE-MODT23DAT\STATUS-NEW’ ) 

OPEN ( UNIT - 14,FILE-’MODT24DAT\STATUS-NEW’ ) 

OPEN (UNIT- 15JTLE-T1JDAT’, STATUS-NEW’ ) 
OPEN(UNIT-16,FTLE-’T2DAT’, STATUS-NEW’) 

OPEN (UNIT- 1 7JFILE«’T3DAT’, STATUS-NEW’ ) 

OPEN (UNIT- 18 JTLE-T4DAT’, STATUS-NEW’ ) 

OPEN (UNIT- 19 JTLE-TSOAT’.STATUS-’NEW’) 
OPEN(UNIT-20TILE-T6X>AT’, STATUS-NEW’) 

OPEN (UNIT-21 .FILE-Nl DAT’, STATUS-NEW ) 

OPEN ( UNIT-22 JTLE-N2DAT’, STATUS-NEW’ ) 

OPEN ( UNIT-23EILE-N30Ar .STATUS-NEW ) 

OPEN ( UNIT -24,FILE-N4DAT’,STATUS-NEW ) 

OPEN (UNIT-25 JP1LE-NSOAT’, STATUS-NEW ) 

OPEN (UNIT-26 JHLE-N60AT, STATUS-NEW ) 


PI-3.1416 

WRITE( 13,* ) EMISSION FOR FOOTPOINTS AND APEX’ 


C*********************************************************** 

C CALCULATE AND WRITE THE KERNAL THERMAL EMISSION AT TIME-0 

C*********************************************************** 
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XALPHA-3.38 

R-1.496E13 

XKAPPA-1.106E-24 

XK-8.617E-8 

XME-9.1E-28 

XMP-1.67E-24 


XN-1.0E11 

T0-2.0E8 

XLO-1.0E8 

AREA-4.132E15 

XI IS THE INTENSITY AT THE EARTH IN PHOTONS PER CM ~ 2 PER KEV 
XEP IS THE ENERGY IN KEV 

CALCULATE XI IN XEP STEPS OF 10 FROM 10 - 100 KEV 


XEP=0.0 

DON=l,21 

XEP-XEP+5.0 

XII- 3.0/((2*PI)**(1.5)) 

C TYPE *,’XI1=’,XI1 

XI2-XKAPPA 
C TYPE*,’XI2=\XI2 

XI3-((XMP/XME)**(.5)) 

C TYPE VXI3*=’,XI3 

XI4=(XN**2)*(AREA/R**2)*(XLO**2)/((XK*TO)**2) 

C TYPE VXI4=’,XI4 

XI5-( (XK*T0/XEP)**4) 

C TYPE *,’XI5-’^Q5 

XI6-EXP ( -XEP / ( XK*T0 ) ) 

C TYPE*,’XI6-\XI6 

XI7=2.0+(2*XEP/ (XK*T0) ) + (XEP**2) / ( (XK*T0)**2) 

C TYPE *,’XI7=’,XI7 

XI8=XI1*XI2*XI3*XI4*XI5*XI6*XI7 
C TYPE *,’XI8-’,XI8 

XI9-XI1*XI2*XI3*XI4*XI5 
C TYPE *,’XI9«\XI9 

XI10«XI9*((.333)*EXP(-XALPHA)*(XEP/<XK*T0))**3) 

C type vxno-’xno 

XII I - XI9* (EXP ( -XALPHA) ) *XI7 

C TYPE*, ’XII 1 -’XH 1 

XI12-XI9* (( -333 ) *EXP( -XALPHA) * (XALPHA**3) ) 

C TYPE *,’XI12-’XI12 

XI(N)-XI8+XI10-XI11-XI12 
C TYPE *,’XI-’XI(N) 

X 1 -XALPHA*XK*T0 
IF(XEP.GTX1)THEN 
XI(N)=0.0 


T83 


C-t> 


ENDIF 

C WRITE( 11,*)’XI-’,XI(N),’KEV«’,XEP 
ENDDO 
C 

C IMS GIVES THE EMISSION FOR A UNIFORMLY TRUNCATED 
C MAXWELLIAN ONLY. TO BE MORE ACCURATE LETS 
C MULTIPLY THIS BY 1 / 3 AND CALCULATE THE EMISSION 
C ASSOCIATED WITH THE DISTRIBUTION PERPENDICULAR TO 
C THIS ONE. WE WILL WEIGHT IT 2/ 3. 

C 

XEP-0.0 

DON1-U1 

XEP-XEP+5.0 

ZIU3.0/((2*PI)**(15)) 

C TYPE* , ZI1-\ZI1 
232-XKAPPA 
C TYPE *,’ZJ2-’JZJ2 

ZI3- ( (XMP/XME) ** ( .5 ) ) 

C TYPE *,’ZI3- , ,ZI3 

234- ( XN**2 ) * ( AREA/R**2 ) * (XLO**2 ) / ( (XK*T0) **2) 

C TYPE *,’ZI4-’,ZI4 

Z35-( (XK*T0/XEP)**4) 

C TYPE ♦.’ZK- ,235 

ZI6-EXP(-XEP/ (XK*TO) ) 

C TYPE *,’ZI6- , ZI6 

ZI7-2.0+(2*XEP/(XK*T0) )+(XEP**2) / ( (XK*T0)**2) 

C TYPE *,’ZT7-’,Zr7 

ZI(N1 )-ZIl*ZI2*ZI3*ZI4*ZI5*ZI6*ZI7 
C TYPE*,’Z3-\ZI(N1),’KEV-\XEP 
ENDDO 
C 

C LET XIT BE TOTAL THERMAL EMISSION IN PH/CM”2/KEV/SEC 
C 

XEP-0.0 

C WRTTE( 13,*) TOTAL THERMAL EMISSION ATTTME-O’ 

C WRircaS.^THOTONS/CM^/KEV/SEC’ 

DOK-1,21 

XEP-XEP+5.0 

XIT(K)«(.333)*XI(K)+(.667)*23(K) 

C TYPE*,’XIT-’XIT(K),’XI- , XI(K), , ZI- , 2I(K) 

C WRITE(13*)’XI-\XI(K), , ZI=\ZI(K) 

XrT(K)-XIT(K)*1.2 

C WRITE(13,*)’KEV-’XEP,Xn'=’XIT(K) 

ENDDO 

C 

£******************«************«******************* 

£*************************************************** 

C CALCULATE EMISSION FROM FOOIPOINTS AT TIME-0 

C*************************************************** 

c*************************************************** 

c 
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CALCULATE GAMMA FUNCTION IN SEPARATE ROUTINE 
REMEMBER TO LINK BOTH MODT2 AND GAMMCP2 TOGETHER 

YK=2.6E- 1 8 !SOME QUESTION HERE 

WRITE ( 13,* ) NONTHERMAL EMISSION ATTTME«=0’ 

XEP=0.0 
DO M= 1,21 
XEP-XEP+5.0 
IF(XEPL,EJC1)THEN 
YI1-(.5)*(PI**(-L5)) 

C TYPE *,’YI1=’,YI1 
YI2=XKAPPA / YK 
C TYPE *, r YI2=’,YI2 

YI3=XN* (AREA/ (R**2) ) *XLO 
C TYPE *, r YI3=’,YI3 

C TYPE *,’XN=’,XN 

C TYPE *,’AREA=’,AREA 

C TYPE *,’XLO=VXLO 

YI4=XK*T0/XEP 
C TYPE *,’YI4=\YI4 

YI5- (XALPHA** ( 1.5 ) ) *EXP( -XALPHA) 

C TYPE *,’YI5=\YI5 

YI6-( ( 1^*XK*T0/XEP)+1 )*GAMMP(XALPHA,15) 

C TYPE *,’YI6=’,YI6 

YI7-YI1*YI2*YI3 
C TYPE*,’YI7= , ,YI7 

YI8=YI7*YI4*YI5 
C TYPE *, r YI8=\YI8 

YI9=YI7*YI6 
C TYPE */YI9=\YI9 

YI(M)=YI8+YI9 

C TYPE *,’XEP=’,XEP,’YI=’,YI(M) 

ENDIF 

IF(XEP.GTJa )THEN 
YI1=(5)*(PI**(-125)) 

C TYPE * ’YIU’.YIl 

YI2 =XKAPP A / YK 
C TYPE *,’YI2-’,YI2 

YI3«XN*AREA*XLO/R**2 
C TYPE*,’YI3-’,YI3 

YI4-((XEP/(XK*T0))**.5)*EXP(-XEP/(XK*T0)) 

C TYPE *, r YI4«’,YI4 

YI5- ( ( 15*XK*T0/XEP) +1.0)*GAMMP(XALPHA,15) 

C TYPE *, r YI5-’,YI5 

YI6-YI 1 * YI2*YI3 
C TYPE *, r YI6«\YI6 

YI7-YI6*YI4 
C TYPE *, r YI7=’,YI7 

YI8-YI6*YI5 
C TYPE *,TI8=’,YI8 
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YI(M)=YI7+YI8 

C TYPE*.’XEP-’,XEP,’YI-’,YI(M) 

ENDEF 

YI(M)-YI(M) *.333 !1 /3 OF TOTAL FLUX COMES IMS WAY 
YI(M)-YI(M)-.3 

C WRITE ( 13,*)’KEV-’,XEP,’YI-’,YI(M) 

C TYPE *.’XEP-’XEP,YI=’,YI(M) 

ENDDO 

C 

C CALCULATE TOTAL EMISSION FROM FLARE AT TIME-0 
C 

XEP-0.0 

C WRITE( 13* ) TOTAL EMISSION AT TIME-0’ 

DOJ-U1 

XEP-XEP+5.0 

TOTI-XIT(J)+YI(J) 

C WRITE( 13,* ) ’KEV-’XEP,TOTI=’,TOTI 

ENDDO 
C 
C 
C 

C******************************************************** 

C CALCULATE EMISSION FROM THERMAL KERNAL AT 1 SEC 
( 2******************************************************** 
C 
C 

C TRUNCATED MAXWELLIAN IN ALONG PARALLEL AXIS 
C WI IS PH/CM' 2/KEV/SEC TRUNCATED 
C RIIS PH/CM'2/KEV/SEC NONTRUNCATED 
C 

XEP-0.0 

XALPHA-2.37 REDEFINED HERE 
T0-6.906E7 TEMP AT 1 SEC 
DOJ-U1 
XEP-XEP+5.0 
WIl-3.0/( (2*PI)**(1.5)) 

C TYPE ‘.’WIl-’.WIl 
WI2-XKAPPA 
C TYPE *,*WI2-’,WI2 

WI3-( (XMP/XME)**(5) ) 

C TYPE *,WI3-\WI3 

WI4-(XN**2)*(AREA/R**2)*(XLO**2) /( (XK*T0)**2) 

C TYPE *,’WI4-’,WI4 

WI5-( (XK*T0/XEP)**4) 

C TYPE *,’W15-’,WI5 

WI6-EXP ( -XEP / ( XK*T0 ) ) 

C TYPE *,WI6-\WI6 

WI7-2.0+(2*XEP/ (XK*TO) )+(XEP**2) / ( (XK*T0)**2) 

C TYPE *,’WI7-’,WI7 

WI8-WI1 *WI2*WI3*WI4*WI5*WI6* WI7 
C TYPE *, r WI8-’,WI8 
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WI9=WI1*WI2*WI3*WI4*WI5 
C TYPE *,’WI9=’,WI9 

WI10=WI9*( (,333)*EXP(-XALPHA)*(XEP/(XK*TO) )**3) 
C TYPE *,’WI10=\WI10 

WI1 1-WI9* (EXP( -XALPHA) ) *WI7 
C TYPE*, ’Will -’.Will 

WI12-WI9* ( ( .333 ) *EXP( -XALPHA) * (XALPHA**3) ) 

C TYPE *,’WI12=\WI12 

WI(J)=WI8+WI10-WI11-WI12 

E1P=14.1 

IF(XEP.0T£1P)THEN 

WI(J)-0.0 

ENDIF 

C TYPE *, r WI=’,WI(J) 

ENDDO 

C 

C REGION TWO CONTRIBUTION 
C 

XEP-0.0 
FCHAR=2.26 1E9 
F1P=2.73E9 
ECHAR-19.3 
E1P-14.1 

SLOP-(FCHAR-FIP) / (ECHAR-E1P) 

YINT=5* (F1P+FCHAR) -,5*SLOP* (ECHAR+E1P) 

DO 1=1,21 
XEP-XEP+5.0 

TI1=(.667)*( (2/XME)**(.5)) 

C TYPE ♦.TH-’, Til 

TI2=XN*XKAPPA*(AREA/(4*PI*(R**2) ) )*XLO/XEP 
C TYPE *,TI2=’,TI2 

TI3=(SLOP*(ECHAR**( 15) ) )+(2*YINT*(ECHAR**(5) ) ) 
C TYPE *,TI3=\TI3 

IF(XEPLE£1P)THEN 

TI4= ( SLOP* ( E1P* * ( 1 5 ) ) ) + ( 2*Y1NT* ( E1P** ( 5 ) ) ) 

C TYPE *,TI4=’,TI4 

ENDIF 

IF(XEP.GT£1PANDJCEPI J T£CHAR)THEN 
TI4= ( SLOP* (XEP** ( 15 ) ) ) + ( YINT* (XEP** ( 5) ) ) 

C TYPE *,TI4=\n4 

ENDIF 

IF(XEP.GT£CHAR)THEN 

TI4-0 

ENDIF 

TI5=-(2*SLOP)/(YINT**2) 

TI(I)=TI1*TI2*TI5*(TI3-TI4) 

C TYPE *,’TI=’,TI(I) 

ENDDO 

C 

C FULL MAXWELLIAN ALONG PERPENDICULAR AXIS 
C AT TIME- 1 SEC 
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c 

XEP-0.0 

DON1-U1 

XEP-XEP+5.0 

RIl-3.0/((2*PI)**(U)) 

C TYPE ♦.’RIU’.RIl 
RI2-XKAPPA 
C TYPE *,’RI2«’,RI2 

RI3-( (XMP/XME)**( 5 ) ) 

C TYPE *,’RI3«’ JU3 

RI4-(XN**2)*(AREA/R**2)*(XLO**2)/((XK*TO)**2) 

C TYPE *.’RI4-\RI4 

RI5-((XK*T0/XEP)**4) 

C TYPE ♦.’RI5-.RI5 

RI6«EXP(-XEP/(XK*T0) ) 

C TYPE *,1116- \RI6 

RI7=2.0+ ( 2*XEP / (XK*T0 ) ) + (XEP**2 ) / ( (XK*TO) **2) 

C TYPE *,’RI7-’JU7 

RI(N1 ) -RI 1 *RI2*RI3*RI4*RI5*RI6*RI7 
C TYPE *,RI=’,RI(N1) ) , KEV=’,XEP 
ENDDO 

LET WIT BE TOTAL THERMAL EMISSION IN PH/CM*2/KEV/SEC 
AT TIME- 1 SEC 


XEP-0.0 

C WRITE (13,*) TOTAL THERMAL INTENSITY AT TIME- 1 ’ 

C WRTTE( 13,*)PHOTONS/CM A 2/KEV/SEC’ 

DOK-U1 

XEP-XEP+5.0 

IF(XEPEE£1P)THEN 

WIT(K)«(.333)*(WI(K)+TI(K))+(.667)*RI(K) 

ENDIF 

IF(XEP.GTEIPAND.XEPLT£CHAR)THEN 

WIT(K)-(.333)*TI(K)+(.667)*RI(K) 

ENDIF 

IF ( XEP.GTECHAR ) THEN 
WIT(K)-.667*RI(K) 

ENDIF 

C TYPE *,’Wrr»\WIT(K),’WI-\WI(K) > RI-’RI(K) 

C WRITE( 13,*)’WI-’,WI(K),’RI«\RI(K) 

Wrr(K)-WIT(K)*lJ 

C WRITE(13,*)’KEV-’ y XEP > r WIT-’ > WIT(K) 

ENDDO 

C 

C********************************************************** 

C CALCULATE NONTHERMAL EMISSION AT TIME- 1 SEC 
£********«************************************************* 
c 

XALPHA-1.12 REDEFINED HERE 
T0-2.0E8 REDEFINED HERE AGAIN 
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XI =XALPHA*XK*TO! REDEFINED HERE AGAIN 
C WRITE (1 3,* iNONTHERMAL EMISSION AT TIME- 1’ 

XEP-0.0 
DOM- 1,21 
XEP-XEP+5.0 
IF(XEPLE.X1)THEN 
SI1=(-5)*(PI**(— 1-5) ) 

C TYPE *,’SI1=’,SI1 

SI2=XKAPPA/YK 
C TYPE *,’SI2«’,SI2 

513- XN* (AREA/ (R**2) ) *XLO 

C TYPE*,’SI3=’,SI3 

C TYPE *,’XN=\XN 

C TYPE *, AREA-’ AREA 

C TYPE *,'XLO=’,XLO 

514- XK*T0/XEP 

C TYPE *,’SI4=\SI4 

515- (XALPHA** ( 15) ) *EXP( -XALPHA) 

C TYPE *,’815=’, SIS 

SI6=( ( 15*XK*T0/XEP)+1 )*GAMMP(XALPHA,15) 

C TYPE *,’SI6=’,SI6 

SI7«SI1*SI2*SI3 
C TYPE *,’SI7=’,SI7 

SI8=SI7*SI4*SI5 
C TYPE *,’SI8«\SI8 

SI9«SI7*SI6 

C TYPE *,’SI9«\SI9 

SI(M)-SI8+SI9 

C TYPE *,’XEP=’XEP,’SI«’,SI(M) 

ENDIF 

IF(XEP.GT.X1)THEN 

SI1=(5)*(PI**(-1.5)) 

C TYPE*,’SI1«\SI1 

512- XKAPPA/YK 

C TYPE *,’SI2=\SI2 

513- XN*AREA*XLO/R**2 

C TYPE *,’SI3-’,SI3 

514- ( (XEP/ (XK*T0 ) ) **5 )*EXP( -XEP/(XK*TO) ) 

C TYPE *,’SI4-\SI4 

515- ( ( 1.5*XK*T0/XEP)+1.0)*GAMMP(XALPHA,15) 

C TYPE *,SI5-\SI5 

516- SI1*SI2*SI3 

C TYPE *,’SI6=’,SI6 

517- SI6*SI4 

C TYPE *,’SI7-’,SI7 

518- SI6*SI5 

C TYPE *,’SI8-’,SI8 

SI(M)-SI7+SI8 

C TYPE *,’XEP«’XEP,’SI«’,SI(M) 

ENDIF 

SI (M ) -SI (M ) *.333 !1 / 3 OF TOTAL FLUX COMES THIS WAY 
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SIX(M)-SI(M)-YI(M) 

SIX(M)-SIX(M)-.3 
C TYPE* , SI-\SI(M),’YI-’,YI(M) 

C WRriE(13,*) , KEV-\XEP,’SK-\SIX(M) 

C TYPE *,’SIX-\SIX(M) 

ENDDO 

C 

C WRITE TOTAL EMISSION AT TIME- 1SEC 
C 

C WRITE (13,*) TOTAL EMISSION AT TIME- 1SEC’ 

C 

XEP-0.0 

DON-1,21 

XEP-XEP+5.0 

TOTIl-SIX(N)+WIT(N) 

C WRITE(13,*)’I£EV»’,XEP,’TOTIl«\TOTIl 
ENDDO 
C 

C********************************************************* 

C CALCULATE EMISSION FOR TIME-2 THRU 30 

C********************************************************* 

c 

C OVERALL LOOP 
C 

XALPHA-3.38 REDEFINED AGAIN 

XLO-I.OE8 

TIME- 1.0 

XKE-1.38E-16 ERO/K 
TEI-2.0E8 [INITIAL TEMP 
DO KT-2,30 

C TYPE * , ’INSIDE LAST LOOP’ 

TIME-TIME+ 1 .0 
C TYPE *,T1ME-\TIME 
TEl-XLO*TEI 

TE2-3* ( (XLO*XKE*TEI/XMP) ** ( 5 ) ) *TIME 
TE3-(XLO)**(15) 

TE-TE1*( (TE2+TE3)**(-.667) ) 

TO-TE 

C TYPE *,TEMP-\T0 
C 

£*****«*************************************************** 
C CALCULATE THERMAL EMISSION FOR TIME-2-30 

C********************************************************* 

C 

C WRITE (13,*) ’APEX THERMAL EMISSION FOR TTME=\TIME 

XEP-0.0 
DON- 1,21 
XEP-XEP+5.0 
VIl-3.0/((2*PI)**(15)) 

C TYPE WIl-’.VIl 

VI2-XKAPPA 
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TYPE *,’V12=’,VI2 
VI3=((XMP/XME)**(.5)) 

C TYPE *, r VI3=\VI3 

VI4«(XN**2)*(AREA/R**2)*(XLO**2)/((XK*TO)**2) 

C TYPE *,’VI4«*\VI4 

VI5«( (XK*T0/XEP)**4) 

C TYPE *, r VI5=’,VI5 

VI6=EXP ( -XEP / ( XK*TO ) ) 

C TYPE *,’VI6=\VI6 

VI7=2.0+(2*XEP/ (XK*TO) )+(XEP**2) / ( (XK*T0)**2) 

C TYPE *,’VI7»’,VI7 

VI8=VIl*VI2*VI3*VI4*VI5*VI6*Vr7 
C TYPE *,’VI8=\VI8 

VI9-VI1*VI2*VI3*VI4*VI5 
C TYPE *,’VI9=’ > VI9 

VI10-VI9* ( ( .333 ) *EXP( -XALPHA) * (XEP/ (XK*TO) ) **3) 
C TYPE *, r VI10-’,VI10 

VII 1-VI9* (EXP ( -XALPHA) )*VI7 
C TYPE *,’VI1 1=\VI1 1 

VI12=VI9* ( ( .333 ) * EXP (-XALPHA) * ( XALPHA* *3 ) ) 

C TYPE *,’VH2=’,VI12 

VI(IO\N)=VI8+VI10-VI11-VI12 
C TYPE *,’VI=’,VI(KT > N) 

XI «XALPHA*XK*TO 
IF(XEP.0TJC1)THEN 
VI(KT,N)-0.0 
ENDIF 

C TYPE */VI«\VI(KT,N) 

C WRnE(ll,*)TffiV«\XEP/VI=\VI(KT,N) 

ENDDO 

C 

C THIS GIVES THE EMISSION FOR A UNIFORMLY TRUNCATED 
C MAXWELLIAN ONLY. TO BE MORE ACCURATE LETS 
C MULTIPLY THIS BY 1/3 AND CALCULATE THE EMISSION 
C ASSOCIATED WITH THE DISTRIBUTION PERPENDICULAR TO 
C THIS ONE. WE WILL WEIGHT IT 2/ 3. 

C 

XEP-0.0 

DON1-U1 

XEP-XEP+5.0 

UIl«3.0/((2*PI)**(l-5)) 

C TYPE *,’UI1«\UI1 

UI2=XKAPPA 
C TYPE *,’UI2«’,UI2 

UI3- ( ( XMP /XME ) ** ( 5 ) ) 

C TYPE *,’UI3-’,UI3 

UI4-(XN**2)*(AREA/R**2)*(XLO**2)/((XK*TO)**2) 
C TYPE *,’UI4=’,UI4 

UI5«((XK*T0/XEP)**4) 

C TYPE *,’UI5=\UI5 

UI6-EXP(-XEP/(XK*T0) ) 
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TYPE *,’UI6«’,UI6 

UI7-2.0+ (2*XEP/ (XK*TO) ) +(XEP**2) / ( (XK*TO) **2) 

C TYPE *,TJT7=’,UI7 

UI(KT,N1 )=UI1*UE*UI3*UI4*UI5*UI6*UI7 
C TYPE*,TJI=\UI(KT,N1),’KEV=’,XEP 

ENDDO 
C 

C LET XTT BE TOTAL THERMAL EMISSION IN PH/CM A 2/KEV/SEC 
C 

XEP-0.0 

C WRITE (13,*) ’APEX THERMAL EMISSION AT ITME= ’.TIME 

C WRITE(13,*)’PHOTONS/CM“2/KEV/SEC 
DO K- 1,21 
XEP-XEP+5.0 

AIT(KTX:) - ( .333 ) *VI(KT JC) + ( .667) *UI(KTJC) 

C TYPE *,’AIT-’^IT(KTJC),’VI-’ > VI(KTJC),’UI- , ,UI(KTJK;) 

C WRTIE(13*)’VI-’,VI(KTJO,Tn-’,UI(KTJC) 
AIT(KT,K)-AIT(KT,K)*1.2 
C WRITE(13,*)’KEV«’XEP ) ’AIT-’AIT(KTX) 

ENDDO 

C 

£*******************************************************«* 

C CALCULATE EMISSION FROM FOOTPOINTS FOR TIME=2-30 SEC 

C********************************************************** 

c 

C WRITE( 13,* ) TOTAL NONTHERMAL EMISSION ATTIME=’,TIME 

XEP-0.0 
DOM-1,21 
XEP-XEP+5.0 
IF(XEPI£JC1)THEN 
BI1-(5)*(PI**(-1J)) 

C TYPE *,311 =\BI1 

BI2-XKAPPA/YK 
C TYPE *,BI2-’3I2 

BI3-XN* (AREA/ (R**2 ) ) *XLO 
C TYPE *,3I3-’3I3 

C TYPE *,’XN«\XN 

C TYPE *,’AREA=’,AREA 

C TYPE *,’XLO«’,XLO 

BI4-XK*T0/XEP 
C TYPE *,3I4«’3I4 

BI5-(XALPHA**( 1.5) )*EXP( -XALPHA) 

C TYPE *,3I5-’3I5 

BI6-((1.5*XK:*T0/XEP)+1)*GAMMP(XALPHA.1.5) 

C TYPE *,3I6-’3I6 

BI7-BI 1 *BI2*BI3 
C TYPE *,317-317 

BI8-BI7*BI4*BI5 
C TYPE *,3I8-’3I8 

BI9-BI7*BI6 

C TYPE *.3I9-’3I9 
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BI(KT,M)=BI8+BI9 

C TYPE *,’XEP=’,XEP,3I=’3I(KT,M) 

ENDEF 

IF(XEP.GTJC1)THEN 
BIl=(3)*(PI**(-l-5) ) 

C TYPE *,3I1=’3I1 

BI2-XKAPPA/YK 
C TYPE *,3I2«’3I2 

BI3-XN*AREA*XLO/R**2 
C TYPE *,3I3=’3I3 

BI4=( (XEP/(XK*T0))**.5)*EXP(-XEP/(XK*T0)) 

C TYPE *,3I4=’3I4 

BI5-( ( 1-5*XK*T0/XEP)+1.0)*GAMMP(XALPHA,13) 

C TYPE *,3I5=’3I5 

BI6=BI 1 *BI2*BI3 
C TYPE *,316= ’316 

BI7«BI6*BI4 

C TYPE *,3I7=’3I7 

BI8-BI6*BI5 

C TYPE *,3I8=’3I8 

BI(KT>1)=BI7+BI8 

C TYPE*.’XEP-’,XEP,3I-’3I(KT,M) 

ENDIF 

BI(KT,M )=BI(KT,M )*.333 !1 /3 OF TOTAL FLUX COMES THIS WAY 
BI(KT,M)=BI(KT,M)-.3 
C WRITE(13,*)’KEV=’XEP,3I=’3I(KT^) 

C TYPE *,’XEP=’,XEP,3I»\BI(KT > M) 

ENDDO 

CALCULATE TOTAL EMISSION FROM FLARE AT TIME=2-30 


XEP=0.0 

WRITE ( 13,*) TOTAL EMISSION AT TIME-’, 
DO J- 1,21 
XEP-XEP+5.0 
TOTI=AIT(KTJ)+BI(KTJ) 

WRITE ( 13,*)’KEV-’,XEP,TOTI=’,TOTI 
ENDDO 

BOTTOM OF OVERALL LOOP 
ENDDO 


********************************************************* 
CALCULATE NONTHERMAL COMPONENT USING THE COLUMN DENSITY 


CALCULATE THE PIXELS WHICH WILL SEE NONTHERMAL EMISSION 
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C CALC LENGTH AND TEMPERATURE ARRAYS 
C 

XKE-1.38E-16 

XLO-1.0E8 

T0-2.0E8 

XMP-1.67E-24 

DO I- 1,21 

nME-(I-l)*1.0 

XXLENl-(3/XLO) 

XXLEN2- (XKE*T0/XMP) ** ( 5 ) 

XXLEN(I)«XLO*( 1+(XXLEN1*XXLEN2*T1ME) )**(.667) 
XXLEN(I)-XXLEN(I)*5 LENGTH 
C TYPE * ' LENGTH OF HALF KERN AL-’^XXLEN(I) 

XXTE(I)«T0*(l+XXLENl*XXLEN2*TIME)**(-.667) .TEMP 
C TYPE*,TEMP«\XXTE(I) 

ENDDO 


C 

C CALC THE LOCATION OF FIRST NONTHERM PIXEL IN TIME 
C 

XNLEN(1)-1.085E9 

XNLEN(2)-1.085E9 

XXPOS(l)-5.422E7 

TOTL-1.1393E9 TOTAL HALF LOOP LENGTH IN CM 
C 

DO JT-2,2I !JT IS TIME IN SECONDS 
DO I- 2,11 

XXPOS(I)-(I-1)*1.085E8 
IF(XXPOS(I)LTJDCLEN(JT) )THEN 
XNLEN( JT) -TOTL-XXPOS(I- 1 ) 

C TYPEVXXPOS(I)«\XXPOS(I) 

ENDIF 

ENDDO 

C TYPE VJT-’JT 

C TYPE*,’XXLEN(JT)«\XXLEN(JT) 

C TYPE VXNLEN(JT)-\XNLEN(JT) 

ENDDO 

C 

C CALCULATE THE COLUMN DENSITY FOR EACH PIXEL IN TIME 
C 

DOKT-1,21 

C 

XFLEN-1.0844E8 ONE LOOP PIXEL LENGTH 
XNCT(KT)-(1.0E12+1.0E13)*XFLEN+ 

1 ( 1.0E11*(XNLEN(KT)-(2*XFLEN) ) ) 

DOKU-1,10 

IF( (KU*XFLEN) EE. (XNLEN(KT)-(2*XFLEN- 1.(^7) ) )THEN 
XNCP (KTKU) -KU*XFLEN* 1 ,0E 1 1 aNDIVPDC COL DENS 
ENDIF 

IF( (KU*XFLEN).GE.(XNLEN(KT)-(XFLEN+1.0E7) ) )THEN 
XNCP(KTJCU)-XFLEN*1.0E13+(XNLEN(KT)-XFLEN)*1.0E11 
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ENDIF 

IF((KU*XFLEN).GE.(XNLEN(KT)-1.0E7))THEN 
XNCP ( KT JOJ ) -XFLEN * ( 1 .OE 1 2 + 1 .OE 1 3 ) + 

1 (XNLEN(KT)-(2*XFLEN))*1.0E11 

ENDIF 

IF( (KU*XFLEN).GE.(XNLEN(KT)+1.0E7) )THEN 
XNCP(KT,KU)«0.0 
ENDIF 

C TYPE *,’KT-\KT,’KU-’jaJ 

C TYPE* PIXEL DENS1TY«’,XNCP(KT,KU) 

C TYPE *,’XNLEN(KT)=\XNLEN(KT) 

ENDDO 

ENDDO 

C 

C 

TYPE*,ENTER PHOTON ENERGY (REAL)’ 

ACCEPT 100.XEP 
100 FORMAT(F8.3) 

TYPE*,ENTERTIME=1,2,3 ) 4,5,6 (INTEGER)’ 
ACCEPT 200 Jvl 
200 FORMAT(I5) 

C 

C PERFORM DOUBLE INTEGRATION 
C CALCULATE AT INFINITY 

C CARVE UP SURFACE FOR SIX TIMES, 0,1,2,5,10,20 
C 

c 

C DOM=l,6 

IF(MEQ.1)THEN 

M7-1 

ENDIF 

IF(M£Q.2)THEN 

M7-2 

ENDIF 

EP(M£Q.3)THEN 

M7-3 

ENDIF 

IF(M£Q.4)THEN 

M7-6 

ENDIF 

IF(M£Q.5)THEN 

M7-11 

ENDIF 

IF(MEQ.6)THEN 

M7-21 

ENDIF 

XNC0L=-.5E19 COLUMN DENSITY INCREMENT 
TYPE *,’XXTE(M7)=\XXTE(M7) 

DO 15-1,1000 
XNCOL-XNCOL+ 1.0E19 
XEO-9.0 !XEO IS THE ELECTRON ENERGY 
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DOJ5-UOO 
XEO-XEO + 2.0 

IF( (XEO**2).OT.(5 JE-18*XNCOL) )THEN 
XJ1«1/(XEP*6*15.0*(6.242E8)**(.5)) 

C TYPE *,’XJ1-’ > XJ1 

C TYPE*,’XEP=\XEP 

XJ2«(XKAPPA*2.0**( U) ) / ( (PI*XME)**(5) ) 

C TYPE *,’XKAPPA-’ ( XKAPPA 

C TYPE *,’XME»’^XME 

C TYPE*,’XJ2-\XJ2 

C TYPE VXK=’ ( XK,’XXTE(M7)=’,XXTE(M7) 

XJ3- ( 1 .0E1 1 ) / ( (XK*XXTE (M7 ) ) ** ( 1 5 ) ) 

C TYPE *,’XK-\XK 

C TYPE *,’XXTE(M7)«\XXTE(M7) 

C TYPE *,’XJ3«\XJ3 

XJ4-XEO*EXP ( -XEO / (XK*XXTE(M7) ) ) 

C TYPE *.’XJ4=\XJ4,’XEO=\XEO 

XJ5-l/(((XEO**2)-(5.2E-18*XNCOL))**(.5)) 

C TYPE *,’XJ5=’ > XJ5 

XI6-1.0E19*2.0 INTEGRATION GRID SIZE 
C TYPE *,’XJ6=’,XJ6 

XJ7-XJ 1 *XJ2*XJ3*XJ4*XJ5*XJ6 
C TYPE *,’XJ7-\XJ7 

XJMX(M7)-XJMX(M7)+XJ7 U AT INFINITY, CLOSE ENOUGH 
C TYPE*,’XJMX(M7)-\XJMX(M7) 

XJEL(I5J5)-XJ7 INTEGRAL DIFFERENTIAL ELEMENT 
ENDS 7 
ENDDO 

C TYPE *,’I5»’45,’J5=M00 

ENDDO 

TYPE *,’M7-'M7 

TYPE *,’XJEL«=’,XJEL(500,500) 


XXEC IS CUTOFF ENERGY 


XXEP IS PHOTON ENERGY PLUS COLUMN DEPTH ENERGY LOSS 

XXEC(M7)-3.38*(3*XK*XXTE(M7) ) 

C TYPEVXXEC(M7)«\XXEC(M7) 

DO 17-1,10 

XXEP(M7) - ( (XEP**2 ) + ( 5 .2E- 1 8*XNCP (M7J7) ) )** ( . 5 ) 

C TYPE*,’XXEP(M7)-\XXEP(M7) 

IF(XXEP(M7).GTJCXEC(M7) )THEN 
XLLIM-XXEP(M7) ! LOWER LIMIT OF INTEGRATION 
ELSE 

XLLIM-XXEC ( M7 ) ! LOWER LIMIT OF INTEGRATION 
ENDS 7 

C TYPE *,’XLLIM«’,XLLIM 
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C TYPE*,’XNCP(M7,I7)=’^NCP(M7 ( I7) 
IF(XNCP(M7J7).GT.0.0)THEN 
XNCOL--.5E19 (COLUMN DENSITY INCREMENT 
DO 18=1,1000 
XNCOL-XNCOL+ 1.0E 1 9 
C TYPE *,’XNCOL=’ v XNCOL 

EF(XNCOLLE JCNCP(M747) )THEN 
XEO=9.0 
DOJ8=UOO 
XEO-XEO+2.0 
IF(XEO.GE JXXLIM )THEN 
XI(M7J7)=XJ(M7J7)+XJEL(I8J8) 

TYPE *,’M7=’M7,’I7= , J7 
TYPE *,’18=’ J8,’J8=’J8 
TYPE *,’XJ(M747) =’ JCJ(M7J7) 

ENDEF 
ENDDO 
ENDIF 
ENDDO 
ENDIF 

TYPE *,’M7=’,M7,T7=\r7 
TYPE *,’XJ(M7 t r7)=\XJ(M7,T7) 

ENDDO 


NORMALIZE J 


DON3=1,10 

IF(XJMX(M7).GT.0.0)THEN 
XXJ(M7,N3)=XJ(M7,N3)/XJMX(M7) NORMJ 
C TYPE *,’XJ (M7.N3 ) =\XJ (M7.N3 ) 

XXN(M7,N3)=5.2E-18*XNCP(M7,N3)/(XEP**2) DIMEN COL DENS 

EF(M7£Q.1)THEN 

WRTIE(21,*)XXN(M7N3) 

WRTIE(22,*)XXJ(M7N3) 

ENDIF 

IF(M7£Q.2)IHEN 
WRITE(23 ,* )XXN(M7N3) 

WRITE(24*)XXJ(M7N3) 

ENDIF 

IF(M7£Q.3)THEN 
WRITE(25 * )XXN(M7N3 ) 

WRTTE(26*)XXJ(M7N3) 

ENDIF 

ENDIF 

ENDDO 


C 

CC CALCULATE I FROM J(I+U+1)-J(U) 
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c 

DOJ5-1.10 

IF(XNCP(M7J5).GT.0.0)THEN 
XXI(M7,1 )«XJ(M7,1 )*4.128E15/ (4*PI* (R**2) ) 

IF(J5.GT.1)THEN 

XXI(M7J5)-ABS(XJ(M7J5)-XJ(M7J5-1 ) ) 

C TYPE*,’M7-’,M7, , J5- , 45 

C TYPE *,’XXI(M7J5)-’4CXI(M745) 

C TYPE *,’XJ(M7,J5)»*’,XJ(M7,J5) 

C TYPEVXJ((M7J5-1)-\XJ(M7J5-1) 

C 

C INTENSITY AT THE EARTH FOR A 1 ARC SECOND CROSS SECT AREA 
C 

XXI(M7J5)=XXI(M7,J5)*4.128E15/(4*PI*(R**2)) 

ENDIF 

C TYPE *,’M7=\M7,’J5=’,J5 

C TYPE *,’XXI(M7J5 ) -’,XXI(M7J5 ) 

ENDIF 
ENDDO 
DO J5- 1,10 
IF(M£Q.1)THEN 
WRITE( 16,* )XXI(M7,1 1-J5) 

ENDIF 

IF(M£Q.2)IHEN 
WRITE( 15,* )XXI(M7,1 1-J5) 

ENDIF 

IF(MJEQ.3)THEN 
WRITE( 17,*)XXI(M7,11-J5) 

ENDIF 

IF(M£Q.4)THEN 
WRITE( 18 ,* )XXI(M7,11-J5) 

ENDIF 

IF(M£Q-S)THEN 
WRTIE( 19,* )XXI (M7,l 1-J5 ) 

ENDIF 

IF(M£Q.6)THEN 

WRITE(20,*)XXI(M7,11-J5) 

ENDIF 
ENDDO 
C FORT-O 
C 

IF(M£Q.1)THEN 

M8-XEP/5 

C TYPE *,’M8«’,M8 

WRITE(16,*)WIT(M8) 

ENDIF 

C 

CFORT-1SEC 

C 

IF(M£Q.2)THEN IT-1SEC 
M8-XEP/5 
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WRITE(15,*)XIT(M8) 

ENDIF 

FORT-2 SECOND 

IF(M£Q.3)THEN IT-2SEC 
MS- XEP /5 

WRITE( 17 ,*)AIT(M7,M8) 
ENDIF 

FORT-5 SECOND 

IF(M EQ.4)THEN !T-6 SEC 
M8-XEP/5 

WRnE(18,*)AIT(M7>i8) 

ENDIF 

FOR T - 10SECOND 

IF(M£Q.5 )THEN !T=11 SEC 
M8-XEP/5 

WRnE(19.*)Arr(M7,M8) 

ENDIF 


FOR T-20SECOND 

IF(M £Q.6 )THEN !T=21 SEC 
M8-XEP/5 

WRITE(20*)AIT(M7Jvl8) 
ENDIF 
C 

DOJ5-1.10 

IF(M£Q.1)THEN 

WRTTE(16,*)XXI(M7J5) 

ENDIF 

IF(MEQ.2)THEN 

WRITE(15*)XXI(M7J5) 

ENDIF 

IF(M£Q.3)THEN 

WRTIE(17*)XXI(M745) 

ENDIF 

IF(M£Q.4)THEN 

WRITE(18*)XXI(M7J5) 

ENDIF 

IF(M£QJ)THEN 

WRTTE(19,*)XXI(M7J5) 

ENDIF 

IF(M£Q.6)THEN 

WRITE(20*)XXI(M745) 

ENDIF 


ENDDO 
C ENDDO 
C 
C 

C NOW FOOTPOINTS INTENSITY ARE AT THE TOP 
C 

C** ************************** ***************************** 

C NOW THAT THE INTENSITIES ARE KNOWN WRITE INTO IDL FILES 
C 

C********************************************************* 

c 

c********************************************************* 

C CLOSING 

c********************************************************* 

c 

TYPE *, PRAISE THE LORD, I AM FREE AT LAST’ 

CLOSE(ll) 

CLOSE(12) 

CLOSE(13) 

CLOSE (14) 

CLOSE(IS) 

CLOSE (16) 

CLOSE(17) 

CLOSE(lg) 

CLOSE(19) 

CLOSE (20) 

CLOSE(21) 

CLOSE (22) 

CLOSE(23) 

CLOSE (24) 

CLOSE(25) 

CLOSE (26) 

C 

STOP 

END 
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PROGRAM MODT 


C ************************************************** 

C THIS PROGRAM CALCULATES THE LOWER CUTOFF VELOCITY OF THE 
C MODEL T ESCAPING ELECTRONS AS A FUNCTION OF TIME 
£******************************************************* 


£******#****************************************************** 
C INITIALIZATION 

£*********** *************************************************** 


c 

VLCO-1.48E10 

TO-O.O 

TTME=0.0 

OPEN ( UNIT- 1 1 JFILE- MODT1 DAT’, STATUS- NEW’ ) 

OPEN ( UNIT= 1 2.FILE- MODT2DAT’ .STATUS- NEW’ ) 
OPEN(UNIT=13 > FILE-’MODT3DAT’, STATUS-NEW’) 

OPEN ( UNIT = 14,FILE=’MODT4 DAT’, STATUS -NEW) 

OPEN (UNIT- 1 5 EILE= MODT5 DAT’ .STATUS- NEW’ ) 

OPEN ( UNIT- 16 JTLE- ’MODT6DAT’,STATUS=’NEW ) 

C 

£*************************************************************** 
C LOOP 

0 ****** ************* ******************************************** 

c 

DO WHILE (IDE. 100) 

I-I+l 


Q ******************************************************** 

C CALCULATE THE LENGTH OF THE KERNAL AT WHICH TIME THE 

C LAST PARTICLE IN THE VLCO FLUX STREAM REACHES THE FRONT 

Q ********************************************************* 

c 

A-(-3.857E12)/VLCO 

B=((-3.857E12)*T0)-1.0E12 

XL— A/ 3.0 

YLl-(A**2)/9.0 

YL2-((-A)**3)/27.0 

YL3«(-B/2.0) 

YL4«((A)**3)*(B/27.0) 

YL5=(B**2)/4.0 

YL6-(YL4+YL5)**.5 

YL7-(YL2+YL3+YL6)**(-1.0/3.0) 

YL-(YL1)*(YL7) 

ZL-(YL2+YL3+YL6)**( 1.0/3.0) 

XLEN- ( XL+ YL+ZL ) * *2.0 
WRITE( 11 ,*)’XLEN-’XLEN 



201 


C CALCULATE NEXT TO" FROM XLEN/VLCO 

C************************************************************ 

c 

TO-TO+ (XLEN/VLCO) 

WRITE (12,*) TO-’, TO 
C 

C*********************************************************** 

C CALCULATE AND SAVE TOTAL TIME ELAPSED 
C************************************************************ 

c 

TIME-TIME+TO 
WRITE( 13,*) TIME-’, TIME 
C 

C************************************************************ 
C CALCULATE TEMPERATURE AT THIS TIME 

C************************************************************* 

c 

TEMP«(2E16)*(((3.857E12)*TIME)+1.0EI2)**(-2.0/3.0) 

WRITE( 14,* )TEMP-’,TEMP 
C 

C************************************************************* 
C SET VLCO AS VMAX, CALCULATE NEXT VLCO 

C************************************************************* 

c 

VMAX-VLCO 

FUN4-2.0E30 

FUN3-1.0E30 

C 

DO WHILE (FUN3LTPUN4) 

C TYPE *, ’VLCO-’, VLCO 

C TYPE *, ’VMAX-’, VMAX 

VLCO - VLCO -1.0E6 
C TYPE *, TEMP-’, TEMP 

A2-(1.024E-6)*(TEMP**(-5) ) 

C TYPE *,’A2-’A2 

B2- ( 3.297E- 12 ) * (TEMP** ( - 1 .0) ) 

C TYPE *,’B2-’32 

Fl-(2.067E15)*(TEMP**(3.0/2.0) ) 

C TYPE *,’F1-\F1 

F2-(A2/B2)*(1.0/B2) 

C TYPE *,’F2-’J 7 2 

F3-(A2*VMAX**2.0) / (2.0*B2) 

C TYPE *,73-’^ 3 

F4-EXP(B2*(VMAX**2.0) ) 

C TYPE *,’F4-’JF4 

F5-EXP(B2*(VLCO**2.0) ) 

C TYPE *,’F5-’J J 5 

F6- ( A2*VLCO**2 ) / ( 2.0*B2 ) 

C TYPE *,’F6-’.F6 

F7-(F2+F3)/F4 
C TYPE *,’F7-’F7 
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FUN1=(F1+F7)*F5 
C TYPE * ) ’FUN1«\FUN1 

FUN2-(F6+F2) 

C TYPE *,’FUN2«’,FUN2 

FUN4=FUN3 

FUN3-ABS (FUN1 -FUN2 ) 
C TYPE *,’FUN4=’J 7 UN4 

C TYPE *,TTJN3-’ JUN3 

ENDDO 
C 

FUN7-VLCO+1.0E6 
WRITE ( 15,*)’VLCO=’,FUN7 
WRTTE( 16* ) TUN3-’TUN3 
C 

ENDDO 



CLOSING 

************************************************************ 


CLOSE(ll) 
CLOSE ( 12) 
CLOSE(13) 
CLOSE ( 14) 
CLOSE (15) 
CLOSE (16) 
C 

STOP 

END 
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APPENDIX C 


SMC End-to-end, Photon Counting Design and Simulation 
This appendix contains an example of the code defining the hard x-ray imag- 
ing telescope described in Chapter IV incorporated into the end-to-end s im ulation 
which gives as an output a data stream of Fourier components. This data stream 
is then processed in AIPS to produce a final image. 
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PROGRAM XSTAR 
C 

C MODELS STARRING XRAY IMAGING TELESCOPE 
C 

C IN THIS MODEL WILL USE TWO SLIT WIDTHS AT 30 DEGREE SPACINGS 
C _ 

C SLIT WIDTHS AND SEPARATIONS WILL CORRESPOND TO THE NEW GRID 

C 

INCLUDE ’XSTARJNC/LISr 
C 

C OVERALL LOOP FOR EACH INDIVIDUAL PHOTON 
C 

C SELECT NUMBER OF PHOTONS PER IMAGE (IPHOT) 

C 

INTEGER*4 EPHOT 

TYPE *,’ ENTER NUMBER OF PHOTONS (INTEGER) PER IMAGE’ 

ACCEPT 100, IPHOT 
100 FORMAT(I) 

TYPE ♦.’ENTER TWIST ANGLE (REAL, IN ARCSECONDS)’ 

ACCEPT 200, TWIST 
200 FORMAT(F8.3) 

TYPE *, TWIST (ARCSEC)-’, TWIST 
TWIST-TWIST* ( 1.0/3600.0) ♦ (PI/ 180.0) 

TYPE ♦, TWIST ( RADIANS ) -’.TWIST 
C 

C IPHOT- 100000 

c 

KUO 
KK3-0 
IDUM— 3 
C 

OPEN(UNIT=3 J 7 ILE=’XSTARDAT’,STATUS=’NEW’) 

OPEN (UNIT- 12 JTLE-’XSTARPOSDAT’, STATUS-NEW’ ) 

OPEN (UNIT- lO.FILE-’XSTARPOINT.DAT’.STATUS-’NEW’) 

C 

C PAINT EXTENDED FLARING ARC 
C 

CALL XSTARFLARE 
C 

C PHOTON LOOP 
C 

DO WHILE (IPHOT.GE.l) 

C 

IPHOT-IPHOT- 1 
PTEST-1 
C 

TYPE *, ’****************’, 'IPHOT-’JPHOT 
C WRITE(3,*) ’IPHOT- ’JPHOT 
C 
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Kl-Kl+1 

C 

CALL XSTARGRID 
C 

C TYPE *,MGRID=’,MGRID,’ ’.’NGRID-7 ’.NGRID 
C WRITE ( 3,* ) ’MGRID- XCRID,’ , r NGRID-\NGRID 
C 

C************GRID ALIGNMENT CHECK*********** 

C 

C MGRID-1 

C NGRID-1 

C 

c****************************************** 

C NOW THAT A GRID HAS BEEN SELECTED, THE PHOTON ARRIVES AT A RANDOM L 
C ON THE GRID WITH A RANDOM ANGLE ASSOCIATED WITH A POINT SOURCE 
C 

C RELATE POINT SOURCE TO TELESCOPE AXIS 
C 

CALL XSTARPOINT 
C 

C TYPE *,’MGRID«’,MGRID,’NGRID«’NGRID 
C 

C RELATE PHOTON ARRIVAL ANGLES TO LOC ALIZED GRID COORDINATE SYSTEM 
C 

CALL XSTARLOCAL 
C 

C TYPE *,’MGRID.’^1GRID,’NGRID=’ > NGRID 
C TYPE *,’THETA3«’,THETA3 
C 

C************GRID ALIGNMENT TEST********* 

C 

C PHI3- 1.45444E-4 

C THETA3-0.0 

C 

C********************************************** 

c 

C DETERMINE POINT OF ARRIVAL ON GRID 
C PUT GRID COORD SYSTEM IN CORNER 
C 

CALL XSTARCOORD 
C 

C TYPE *,MGRID= , > MGRID, , NGRID= , ,NGRID 
C 

C *********UN!F OR M PHOTON DISTRIBUTION TEST************* 

C 

C XPOS-(10.0/10000)*(IPHOT+1) 

C YPOS-O.O 

C 

c****************************************************** 

c 

C START SECOND LOOP TO TEST FOR PHOTON THROUGHPUT SURVIVAL 
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TEST FOR PHOTON PENETRATION THROUGH THE FIRST GRID 
CALLXSTAR PEN 


TYPE *, ’MGRID=’,MGRID,’NGRID=’,NGRID 
TYPE *, ’PTEST-’FTEST 


IF (PTEST.EQ.1 ) THEN 
CALL XSTARXY 

TYPE *,MGRID= , MGRID,’NGRID-’ I NGRID 


CALL XSTARPEN2 

TYPE *,’MGRID=’ > MGRID,’NGRID=’ > NGRID 
ENDIF 

TYPE *, ’PTEST-’JPTEST 

DETERMINE PHOTON POINT OF IMPACT ON DETECTOR 
IF (PTESTLQ.l) THEN 
CALL XSTARXY2 

TYPE * .’MGRID-’MGRID.’NGRID- ’.NGRID 
DETERMINE FOURIER COMPONENT, BUILD MATRIX 
CALL XSTARFOUR 
WRITE XPOS MATRIX 

CALL XSTARMTX (KX3) 

WRITE( 12* )XPOS,YPOSXPOS2,YPOS2XPOS3,YPOS3 
ENDIF 

ENDDO 

WRITE FOURIER COMPONENT MATRIX TO AEPS FORMAT 

TYPE *,’NTOT=’NTOT 
WRITE ( 3,* ) ’NTOT-VNTOT 
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DONl-l.TSLIT 
DO N2-1.TGRID 

WRITE(3,* ) MGRID-’,Nl,NGRID-’,N2 > DTOT- , DTOT(Nl,N2) 
WRITE(3,*)BIN(N1,N2,1)3IN(N1,N2,2)3IN(N1,N2,3), 

1 BIN(N1,N2,4)3IN(N1,N24)3IN(N1,N2,6), 

1 BIN(N1,N2,7)3IN(N1,N2,8)3IN(N1,N2,9), 

1 BIN(N1,N2,10)3IN(N1,N2,11)3IN(N1,N2,12). 

1 BIN(N1,N2,13)3IN(N1,N2,14)3IN(N1,N2,15), 

1 BIN(N1.N2,16)3IN(N1,N2,17)3IN(N1,N2,18), 

1 BIN(N1,N2,19)3IN(N1,N2,20)3IN(N1,N2,21) 

ENDDO 

ENDDO 

C 

DON3-1.21 

WRITE(10,*)’N3-’,N3, , PANG(N3)-\PANG(N3) 

WRITE( 10 *)THETA(N3)- , ,THETA(N3) 
WRITE(10,*)’PHI(N3)- , 3HI(N3) 

ENDDO 

C 

CALL XSTARFMAT 
C 

C PLOT PHOTON FLUX VERSUS X ANDY 
C 

C CALL XSTARPLOT 

C 

c 


STOP 

END 

C 

C XSTARJNC 

REAL*4 DIV, NDIV, MDIVRAND 
C 

C 

C 

C 

C 

C 

C 

C 

C 

C 


REAL*4 MDIV2.NDIV2 

REALM DELTHETA(50),NUM,DENOM 

REAL*4 PSI(55 ) ,THETA(55 ) ,PHI(55) 

REAL*4 DELTHL( 1045) . SLITW( 10) ,SLITW2( 10) 
REAL*4 XPOS,YPOS ( XPOS2,YPOS2 t XPOS3,YPOS3 
REAL*4 SHIFT ( 10) .STIF.GRIDW 
REAL*4 SHIFT2( 10) ,SNUM 
REAL*4 VX( 110,110) ,VY( 1 10,1 10) 

REAL*4 XDET( 540,30) 

REAL*4 DTOT( 1045 ) 3EM 
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c 


c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 


REALM U( 110, 110), V( 110,110) 

REALM DANGJDPOSJDETRESJDETANG 
INTEGERS TGRID, MGRIDNGRID,K2 
INTEGERM PTEST, TSLIT 
INTEGERM ISNUM, NSLIT 
INTEGER*4 KSNUM, PANG ( 30 ) 

INTEGERM NTOT.INCR3IN (5,50,30 ) 

INTEGERM REM2 

COMMON /GRID/ MGRIDNGRID,ISEED,SLITW 

COMMON /GRIDC / TGRID, TSLIT 

COMMON /POINT/ THETA3 ,PHI3 PANG ,THETA,PHI 

COMMON/COORD/ XPOS.YPOS 

COMMON /LOC/ DELTHL 

COMMON /PTEST/ PTEST 

COMMON /XY/ XPOS2.YPOS2 

COMMON /PEN2 / ISNUM 

COMMON /XY2 / XPOS3.YPOS3 

COMMON/FOUR/ VX.VY 

COMMON / MTX / BIN ,XDET ,NTOT JDTOT 

DATAPI/3.14159/ 

DATA ISEED/ 1 / 

REALM FUNCTION RANO(IDUM) 

JUST SET IDUM TO ANY NEGATIVE SEED TO START 


DIMENSION V(97) 

DATA IFF /0/ 

IF ( IDUM -LT.0. OR. IFF EIQ.O) THEN 
IFF-1 


ISEED= ABS ( IDUM ) 

IDUM-1 

DO 11 J=l,97 
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DUM-RAN(ISEED) 

11 CONTINUE 
DO 12 J- 1,97 

V(J)-RAN(ISEED) 

12 CONTINUE 
Y-RAN(ISEED) 

ENDIF 

J=1+DMT(97.*Y) 

IF( J.GT.97.0R JJLT. 1 )PAUSE 
Y-V(J) 

RANO-Y 

V(J)-RAN(ISEED) 

RETURN 
END 

SUBROUTINE XSTARCOORD 

THIS ROUTINE DETERMINES THE RANDOM LOCATION ON THE SURFACE OF 
THE FIRST GRID ON WHICH THE PHOTON FALLS 

INCLUDE ’XSTARJNC/LIST’ 

TYPE *, ENTERING XSTAR COORD’ 

RAND-RANO ( IDUM ) 

TYPE ♦.’RAND-’JIAND 

XPOS-RAND* 10 ON CM 

RAND-RANO (IDUM) 

TYPE *, ’RAND- ’.RAND 

YPOS=RAND*10 HNCM 

TYPE*,’XPOS-’XPOS 
TYPE *,YPOS-’,YPOS 

RETURN 
END 

SUBROUTINE XSTARFLARE 
C 

C THE PURPOSE OF THIS SUBROUTINE IS TO BUILD A FLARE 
C FROM A COLLECTION OF POINT SOURCES, A SIGNIFICANT 
C PROBLEM IS WHAT COLLECTION OF POINTS BEST REPRESENT 
C A FLARE OF INTEREST. 

C 

C THIS ATTEMPT AT PAINTING A FLARE WILL BE REPRESENTED 
C BYA 1 00 ARCSECOND LONG SEMICIRCLE WITH 20 POINT 
C SOURCES EQUA- SPACED ALONG THE ARC. THIS MEANS THE 
C RADIUS OF THE FLARE WILL BE 100/PI. 

C 

C THE ORIENTATION OF THE FLARE TO THE VERTICAL AXIS 
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C OF THE TELESCOPE WILL BE 45 DEGREES. THE FIRST 
C POINT WILL BE 0.0 DEGREES OF ROTATION FROM THE 
C HORIZONTAL, THE NEXT WILL BE 9.0, THE NEXT 
C 9.0+9.0 AND SO FORTH. ESSENTIALLY, DIVIDING 
C 180 DEGRESS BY 21 PIXELS AND PUTTING THE POINT 
C SOURCE IN THE CENTER OF EACH PIXEL. 

C 

INCLUDE ’XSTARJNC /LIST’ 


TYPE *, ’ENTERING XSTARFLARE’ 

WRITE ( 10,* ) 'ENTERING XSTAR FLARE’ 

LET PSI BE THE ANGLE AROUND THE CIRCLE 


1-0 

PSI(l) -0.0 

C TYPE*,’PSI(1)-\PSI(1) 

DELPSI-PI/20 
DO I = 2,21 

PSI(I)-PSI(J-1)+(DELPSI) 

C TYPE *,’I-M,’PSI(I)=’,PSI(I) 

END DO 
C 

C DEFINE ANGLES OF FLARE POINT 
C 

C ANGLE OF ROTATION ABOUT Z-AXIS OF TELESCOPE 
C 

C CENTER FLARE AT 45 DEGREES ABOUT CENTER 
C 30 ARCSEC FROM CENTER OF TELESCOPE 
C 

CLOOP 

C 

12-0 

C 

DO WHILE (J2TE.20) 

C 

J2-I2+1 

C 

PHI2- 30 'ARCSECONDS FROM CENTER OF TELESCOPE 
C PHI2-60 (ARCSECONDS FROM CENTER OF TELESCOPE 
C PHI2-240 [ARCSECONDS FROM CENTER OF TELESCOPE 
C PHL2-480 (ARCSECONDS FROM CENTER OF TELESCOPE 
C PHI2-960 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2-540 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2-600 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2-700 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2-800 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2-900 (ARCSECONDS FROM CENTER OF TELESCOPE 

C PHI2- 1000 (ARCSECONDS FROM CENTER OF TELESCOPE 

C 
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C TYPE******************************************* 

C TYPE ♦.’ENTER ARCSECONDS OUT FROM SUN CENTER (INT ) ’ 

C TYPE *,’******************************************’ 

C ACCEPT 1000, PHI2 
C1000 FORMAT(I) 

C 

ARCRAD-10.0 

C 

NUM-(ARCRAD)*COS(PSI(J2)) UN ARCSECONDS 
C TYPE *,’NUM»\NUM 

DENOM- ( (ARCRAD)*(SIN(PSI(J2) ) ) )+PHI2 ON ARCSECONDS 
C TYPE *,DENOM-’DENOM 

DELTHETA ( J2 ) -AT AN ( NUM /DENOM ) UN RADIANS 
C 

THETA(J2)-( (45.0)* (PI/180) )+DELTHETA(J2) 9N RADIANS 
TYPE *,’J2«M2.* \THETA(J2)«\THETA(J2) 

ENDDO 

ANGLE OF ROTATION AWAY FROM Z-AXIS OF TELESCOPE 
LOOP 
J3-0 

TYPE*,’PSI(1)-\PSI(1) 

DO WHILE (J3LE.20) 

J3-J3+1 

TYPE *,’J3-\J3 
TYPE *,’PSI(J3)- , J»SI(J3) 

TYPE *,?HI2-\PHL2 

NUM-PHI2+( (ARCRAD)*SIN(PSI(J3) ) ) 

NUM-NUM*(4.848E-6) CONVERTING PHI TO RADIANS 
DENOM -COS ( DELTHETA ( J3 ) ) 

PHI(J3 ) -NUM /DENOM 
TYPE *,’DELTHETA(J3)-\DELTHETA(J3) 

TYPE *,’NUM-\NUM,’ \denom-\denom 

TYPE *,’J3-V 7PHI(J3)-\PHI(J3) 

ENDDO 

TYPE *,’PHI3- , THI3,’ ’,THETA3-’,THETA3 
WRITE(10*) T>HI3«\PHI3,’ ’.THETAS-’.THETAJ 

RETURN 
END 
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SUBROUTINE XSTARFMAT 
C 

C THIS ROUTINE WRITES UV FILE FOR AIPS /UVFIL 

C 

INCLUDE ’XSTAR JNC /LIST’ 

C 

TYPE *,’ ENTERING XSTAR FMAT’ 

C 

OPEN(UNIT=4 > FILE=’IN2FILX)AT’,STATUS=’NEW’) 

C 

DOJI-l.TSLIT 

DOJ2-1.TORID 

Rl-SLITW(Jl) 

R2=520.0 !GRID SEP IN CM 
R3=R2/R1 

U(J1J2)=R3* (SIN(DELTHL( J1J2) ) ) 

V(JlJ2)-R3*(COS(DELTHL(JU2))) 

C 

U(JU2)=-U(J142)/(1.0E4*2.0) 110-5-90 
C U(J1J2)-U(JU2)/(1.0E4*PI) !OLD 

C 

V(JU2)-+V(JU2)/(1.0E4*2.0) 110-5-90 
C V(J1J2)— V(J1J2)/(1.0E4*PI) !OLD 

C 

VX(JU2)— VX(JU2) LAJPS CONVENTION 
C 

VY(J1J2)— VY(J1J2) 1AIPS CONVENTION 
C 

WT=DTOT(JU2) 

WRnE(4,2001)U(JU2),V(IU2),VX(JU2),VY(JU2).WT 
2001 FORMAT(2XE10.6,1XE10.6^XJ J 12.6,1XE12.6,1XE10.4) 

ENDDO 

ENDDO 

C 

RETURN 

END 

SUBROUTINE XSTARFOUR 
C 

C THIS ROUTINE CALCULATES THE FOURIER COMPONENTS AND WRITES THEM 
C TOAFILE IN THE APS/UVFIL FORMAT 
C 

INCLUDE ’XSTAR JNC /LIST’ 

C 

TYPE *, ENTERING XSTAR FOUR’ 

C 

C DETECTOR SPATIAL RESOLUTION IS IMPORTANT HERE 
C 

c 

C DANG-((XPOS3*2*PI)/10)-(PI) [RADIANS 
DANG- (XPOS3*2*PI) / 10 
C 
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M-MORID 

N-NORID 

VX(M,N)-VX(M > N)+COS(DANG) 

VY(M,N ) -VY(MN ) +SIN(DANG ) 

C 

C TYPE ♦,7^GRE)- , > MGRID,’NGRID-\NGRID 
C TYPE ‘.M-’.M.’N-’N 
C TYPE *,’VX(M,N)«\VX(M,N) 

C TYPE *,’VY(M,N) -\VY(M,N) 

RETURN 

END 

SUBROUTINE XSTARPEN 
C 

C THIS SUBROUTINE DECIDES WHETHER THE PHOTON SURVIVES STRIKING 
C THE FIRST GRID 
C 

INCLUDE XSTARJNC/LIST’ 

C 

TYPE *, ’ENTERING XSTARPEN’ 

C WRITE (3,*) ’ENTERING XSTAR PEN’ 

C 

C SLITW is GIVEN IN XSTAR GRID WHEN GRID IS CHOSEN 
C 

M-MGRID 

C 

REM-AMOD( 10.0,SLITW( 1 ) ) 

SHIFT( 1) -REM/2 
SHEFT(2) -.0125/2 
C TYPE *,’REM-’,REM 
C TYPE *, ’SHIFT-’, SHIFT 
C TYPE VISNUM-’JSNUM 
C 

C CENTER GRIDS IN 10 CM SQUARE AREA 
C 

CLOOP 

C 

C FIND WHICH INTERVAL SLIT OR SLAT THAT THE PHOTON FALLS 
C 

Jl-0 

STIF-0.0 

DO WHILE (XPOS.GT.STIF) 

Jl-Jl+1 

C TYPE *,’J1-’J1 

STIF-(J1*SLITW(M) )+SHIFT(M) 

C TYPE *,’STIF-’,STIF 
NSLIT-J1 
ENDDO 

C TYPE ‘.’NSLIT-’NSLIT 
C 

C DETERMINE IF IT PASSES THROUGH 
C 
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REM2- MOD(NSLIT^) 

C TYPE *,’REM2«=’,REM2 
IF (REM2.GT.0) THEN 
PTEST-7 
ENDIF 
C 

RETURN 

END 

SUBROUTINE XSTARGRID 
C 

C THE PURPOSE OF THE SUBROUTINE IS TO RANDOMLY SELECT A GRID 
C FOR PHOTON PENETRATION 
C 

C IN THIS CASE TWO SLITS WIDTHS AT 30 DEGREES SEPARATION LEADS TO 
C 24 GRIDS 
C 

INCLUDE ’XSTAR JNC /LIST’ 

C 

C WRITE( 3 * ) ENTERING XSTAR GRID’ 

C TYPE * .ENTERING XSTAR GRID’ 

C 

TSLIT-2 

TGRID-24 

C 

RAND=RANO ( EDUM ) 

C TYPE *,’RAND=\RAND 
C WRITE (3,*) ’RAND=’,RAND 
MDIV-1.0/TSLIT 
NDIV= 1 .0 /TGRID 
K2-0 
C 

C 1ST INDEX 
C 

DO WHILE (K2LE.100) 

K2-K2+1 

C TYPE *,’K2«\K2 
MDIV2=K2*MDIV 
C TYPE *,’MDIV2-’,MDIV2 
IF (RANDLEAOIV2) THEN 
MGRID=K2 
K2-101 
ENDIF 
ENDDO 
C 

C 2D INDEX 
C 

K3-0 

RAND-RANO(IDUM) 

C TYPE ♦.’RAND-’JIAND 
DOWHILE (K31E.100) 
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K3-K3+1 
NDIV2-K3*NDIV 
IF ( RAND-LE-NDIV2 ) THEN 
N0RID-K3 
K3-101 
ENDIF 
ENDDO 
C 

C DETERMINE SLIT WIDTH HERE 
C 

SLITW(l)-.0275 

SLITW(2)-.0125 

C 

C TYPE * .’LEAVING XSTAR ORID’ 

C 

RETURN 

END 

SUBROUTINE XSTARLOCAL 
C 

C THIS ROUTINE ROTATES THE INCOMING PHOTON ANGLE TO CORRESPOND TO T 
C OF THE GRID IN THE TELESCOPE 
C 

INCLUDE XSTARJNC/LIST’ 

C 

C THIS ROUTINE MUST ALSO ASSIGN ROTATION ANGLES TO ALL GRIDS 
C 

TYPE *, "ENTERING XSTARLOCAL’ 

C 

DELTHL (1,1) >0.0 
DELTHL(2,1)«0.0 
DOJ1-1.TSLIT 
DOJ2-2.TGRED 

DEL1HL(JU2)-DELTHL(J1 J2-1 )+( 15*(PI/180) ) 

C TYPE*’JU\J1,’7J2-\I2 
C TYPE*,’DEL1HL(JU2)’T)ELTHL(JU2) 

ENDDO 

ENDDO 

C 

C TYPE ♦.’DELTHL’ 

C TYPE*, DELTHL 
C 

C APPLY LOCAL GRID ANGLE TO INCOMING PHOTON ANGLE ABOUT TELESCOPE A 
C 

THETA3=THETA3-DELTHL ( MGRID.NGRID ) 

C 

RETURN 

END 

SUBROUTINE XSTARMTX (KK3 ) 

C 

C THIS ROUTINE GIVES THE TOTAL NUMBER OF PHOTONS MAKING IT 
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C AND THE TOTAL NUMBER MAKING IT TO A DETECTOR 
C AND THE NUMBER PER BIN 
C 

INCLUDE ’XSTAR.INC /LIST’ 

C 

INCR-21 

DETRES-10.0/INCR 

KK3-KK3+1 

C 

Jl-0 

DPOS-O.O 

DOWHILE(XPOS3.GTDPOS) 

Jl-Jl+1 

DPOS-J 1 *DETRES 
ENDDO 
C 

M-MGRID 

N-NGRID 

C TYPE *,’J1=’J1 

BIN(M,NJ1)=BIN(MW1)+1 
XDET ( M.N J 1 ) =DPOS -DETRES / 2 
NTOT-KK3 

DTOT(MJN) -DTOT(M,N) + 1.0 

TYPE *,’NTOT=’J*JTOT 

RETURN 
END 

SUBROUTINE XSTARPEN2 

THIS SUBROUTINE DECIDES WHETHER THE PHOTON SURVIVES STRIKING 
THE FIRST GRID 

INCLUDE ’XSTARJNC /LIST’ 

TYPE *, ’ENTERING XSTARPEN2’ 

WRITE ( 3.* ) ENTERING XSTARPEN2’ 

SLITW IS GIVEN IN XSTAR GRID WHEN GRID IS CHOSEN 
CENTER GRIDS IN 10 CM SQUARE AREA 

M-MGRID 

REM-AMOD( 10.0.SLITW ( 1 ) ) 

SLITW2 ( 1 ) = ( 10.0-REM) / ( 365 ) 

SLTTW2(2)-( 10.0-.0125) / (801) 

SHIFT2(1) -(REM/2) 

SHIFT2( 2) -.0125/2 

NOTE DIFF HERE FROM FIRST ROUTINE 
TYPE *,’REM-\REM 
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C TYPE *,’SLITW2=\SLITW2 
C 

CLOOP 

C 

C FIND WHICH INTERVAL SLIT OR SLAT THAT THE PHOTON FALLS 
C 

11-0 

STIF-0.0 

DO WHILE (XPOS2.OT.STIF) 

Jl-Jl+1 

C TYPEVJl-’Jl 

STIF-(J1*SLITW2(M) )+SH3FT2(M) 

C TYPE *,’SnF=’,STIF 
NSLIT-J1 
ENDDO 
C 

C TYPE ♦.’NSLIT-’.NSLIT 
C 

C DETERMINE IF IT PASSES THROUGH 
C 

REM2=MOD(NSLIT,2) 

C 

C TYPE *,’REM2«’,REM2 
C 

IF (REM2.EQ.0) THEN 
PTEST-7 
ENDBF 
C 

RETURN 

END 

SUBROUTINE XSTARPOINT 
C 

INCLUDE ’XSTARJNC/LIST’ 

C 

c 

TYPE ♦.’ENTERING XSTAR POINT’ 

C TYPE *, ’CONFIRMED’ 

C WRITE( 10 * ) ’ENTERING XSTAR POINT’ 

C 

C NOW RANDOMLY SELECT ANGLES FOR PHOTON 
C 

RAND-RANO(IDUM) 

C 

C FOR UNIFORM FLARE OF 21 POINT SOURCES 
C 

C DIV-1.0/21 
C 
C 
C 

TYPE ♦.’RAND-’RAND 
C 
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c 

CLOOP 

c 

C J4-0 
C NDIV=0.0 

c 

C DOWHILE (RAND.GT.NDIV) 

C 

C J4-J4+1 
C TYPE VJ4=’J4 
C NDIV=J4*DIV 
C TYPE *,’NDIV-’ > NDIV 
C Pffl3=PHI(J4) 

C TYPE *,’PHI3=’,PHI3 

C THETA3-THETA(J4) 

C TYPE *,THETA3=’,THETA3 

C 

C ENDDO 
C 

C FOR UNIFORM FLARE OF 21 POINT SOURCES 
C WITH ENDPOINT 10 TIMES BRIGHTER 
C 

C DIV=1.0/39 
C 

C P6-DIVM0 
C IF (RANDiEF6)THEN 
C PH13-PHI(1) 

C THETA3=THETA( 1 ) 

C 16=1 
C ENDIF 
C 

C P7=DIV*29 
C IF (RAND.GEP7)THEN 
C PHI3=PHI(21) 

C THETA3=THETA(21 ) 

C J6=21 
C ENDIF 
C 

C 12=0 

C IF (RAND.GTP6ANDHANDLrTP7)THH>I 
C DOJS-2,20 
C 12=12+1 

C P1=((DIV*10)+(I2*DIV) ) 

C IF ( RAND. GT P6 AND PANDLEP1 )THEN 

C PHI3=PHI(J5) 

C THETA3=THETA(J5) 

C J6=J5 

C ENDIF 
C P6-P1 
C ENDDO 
C ENDIF 
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cc 

C TYPE *,’J4-\J4 

C 

C 

C PANG(J4)-PANG(J4)+1 
C 

C PANG (J6) -PANG (J6)+l 
C 

C TYPE J5 

C TYPE *,’J6-’ J6.TANG (J6) J»ANG( J6) 

C TYPE* HAND-’ JIAND 
C 

C TYPE *,THI3-’j > HI3,’ ’,THETA3-’,THETA3 
C WRITE ( 10.* ) THI3=’Pffl3.’ ’,THETA3-’,THETA3 
C 

C****************************************************************** 
C***»* FOR FLARE WITH ONLY TWO ENDPOINTS /FOOTPODMTS*** 

£******************************************4 *************^ ********* 

c 

C TYPE *,’***********************************• 

C TYPE *,TWO FOOTPOINT FLARE’ 

C TYPE *’*****»*****************************• 

C 

C IF (RAND.LT..5)THEN 
C PHI3-PHI(1) 

C THETA3-THETA( 1 ) 

C J8-1 

C ENDIF 

C IF (RAND.GE.J)THEN 
C Pffl3=PHI(21) 

C THETA3-THETA(21) 

C J8-21 

C ENDIF 

C PANO(J8)-PANG(J8)+l 
CC 

£******************************************************* 

C*** FOR FLARE WITH WORST CASE PHOTON VARIATION******** 

£**«**************************************************** 

c 

C XNA-4900.0 
C XNB -4624.0 
C XNTOT-99864.0 
C XNRATA-XNA / XNTOT 
C TYPE *,’XNRATA-’,XNRATA 
C XNRATB-XNB /XNTOT 
C TYPE *,’XNRATB-\XNRATB 
C IF (RANDiTJCNRATAANDJlAND.GE.0.0) THEN 
C PHI3-PHI(1) 

C THETA3 -THETA ( 1 ) 

C PANG(1)-PANG(1)+I 

C TYPE*RANG(1) 
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C ENDIF 

C IF (RAND1T.(2*XNRATA) AND.RAND.GE.XNRATA) THEN 
C PHI3=PHI(2) 

C THETA3<=THETA(2) 

C PANG(2)-PANG(2)+1 

C TYPE *,PANG(2) 

C ENDIF 

C IF (RAND.LT.(3*XNRATA) ANDJtAND.GE.(2*XNRATA)) THEN 

C PHI3=PHI(3) 

C THETA3=THETA( 3 ) 

C PANG(3)=PANG(3)+1 

C TYPE*JPANG(3) 

C ENDIF 

C IF (RAND1T.(4*XNRATA)AND.RAND.GE.(3*XNRATA) ) THEN 
C PHI3-PHI(4) 

C THETA3=THETA(4) 

C PANG(4)-PANG(4)+1 

C TYPE *,PANG(4) 

C ENDIF 

C IF (RAND.LT. (5*XNRATA ) ANDHAND.GE. (4*XNRATA) ) THEN 
C PHI3=PHI(5) 

C THETA3-THETA(5) 

C PANG(5)-PANG(5)+1 

C TYPE *,PANG ( 5 ) 

C ENDIF 

C XNR1=(5*XNRATA) 

C TYPE *,’XNR1=\XNR1 
C XNR2=XNR1 +XNRATB 
C TYPE *,’XNR2«’,XNR2 

C IF ( RAND JLTJCNR2 AND JIAND.GEJCNR1 ) THEN 
C PHI3=PHI(6) 

C THETA3=THETA(6) 

C PANG(6)=PANG(6)+1 

C TYPE*PANG(6) 

C ENDIF 

C XNR3-XNR1+(2*XNRATB) 

C TYPE *,’XNR3-’,XNR3 

C IF ( RAND LT.XNR3 AND .RAND.GE.XNR2) THEN 
C PHI3-PHI(7) 

C THETA3=THETA ( 7 ) 

C PANG(7)-PANG(7)+1 

C TYPE *J*ANG(7) 

C ENDIF 

C XNR4=XNR1+(3*XNRATB) 

C TYPE *,’XNR4-’,XNR4 

C IF (RANDLT-XNR4AND.RAND.GE.XNR3 ) THEN 
C PHI3=PHI(8) 

C THETA3»THETA(8) 

C PANG(8)-PANG(8)+1 

C TYPE * J*ANG ( 8 ) 

C ENDIF 
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C XNR5«XNR1+(4*XNRATB) 

C TYPE *,’XNR5-\XNR5 

C IF (RANDJJT.XNR5AND.RAND.GE.XNR4) THEN 
C PHI3«PHI(9) 

C THETA3-THETA( 9 ) 

CC PANG(9)-PANG(9)+1 
CC TYPE*J»ANG(9) 

C ENDIF 

C XNR6«XNR1+(5*XNRATB) 

C TYPE *,’XNR6-’,XNR6 

C IF ( RANDiT JCNR6 AND JIAND.GE JG4R5 ) THEN 
C PHI3-PHI(10) 

C THETA3 -THETA ( 10 ) 

C PANG(10)«PANG(10)+1 

C TYPE *J*ANG( 10) 

C ENDBP 

C XNR7-XNR1+(6*XNRATB) 

C TYPE VXNR7-.XNR7 

C IF (RANDJLTJCNR7AND.RAND.GE.XNR6) THEN 
C PHI3-PHI(11) 

C THETA3 -THETA ( 1 1 ) 

C PANG(11)-PANG(11)+1 

C TYPE * J*ANG( 1 1 ) 

C ENDIF 

C XNR8-XNR1+(7*XNRATB) 

C TYPE *,’XNR8-’,XNR8 

C IF ( RAND iT_XNR8 AND JRAND.GEJXNR7) THEN 
C PM3«PHI(12) 

C THETA3-THETA( 12) 

C PANG(12)=PANG(12)+1 

C TYPE *,PANG( 12) 

C ENDIF 

C XNR9-XNR1+(8*XNRATB) 

C TYPE *,’XNR9-’,XNR9 

C IF ( RAND.LT JCNR9 AND JIAND.GE JCNR8 ) THEN 
C PHI3-PHI(13) 

C THET A3 - THETA ( 13 ) 

C PANG(13)-PANO(13)+l 

CC TYPE*PANG(13) 

CC ENDIF 

CC XNR10-XNR1+(9*XNRATB) 

CC TYPE *,’XNR10-’,XNR10 

CCC IF ( RAND I.T JXNR 1 0 AND JIAND.GEJCNR9) THEN 
C PHI 3 -PHI ( 14) 

C THETA3-THETA ( 14 ) 

C PANG(14)-PANG(14)+1 

C TYPE*,PANO(14) 

C ENDIF 

C XNR11-XNR1+(10*XNRATB) 

C TYPE*/XNRll-’ r XNRll 

C IF ( RAND.LT.XNR 1 1 AND.RAND.GE-XNR 1 0 ) THEN 
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C PHI3=PHI(15) 

C THETA3-THETA( 15 ) 

CC PANO(15)-PANG(15)+l 

C TYPE *J > ANG( 15) 

CC ENDIF 

C XNR12-XNR1+(11*XNRATB) 

C TYPE *,’XNR12=’^XNR12 

C IF (RAND.LT JCNR12 AND RAND.GE.XNR1 1 ) THEN 
C PHI3-PHI(16) 

C C THETA3*=THETA( 16) 

C PANG( 16)=PANG( 16)+1 

CC TYPE*J > ANG(16) 

C ENDIF 

C XNR13-XNR12+XNRATA 
CC TYPE *,’XNR13-’,XNR13 

CC IF (RANDTT_XNR13ANDRAND.GEJXNR12) THEN 
C PHI3=PHI( 17) 

C THETA3-THETA(17) 

C PANG(17)-PANG(17)+1 

C TYPE *,PANG( 17) 

C ENDIF 

C XNR14«XNR12+(2*XNRATA) 

C TYPE *,’XNR14»’,XNR14 

C IF (RANDTTJCNR14ANDRAND.GEJCNR13) THEN 
C PHI3-PHH18) 

C THETA3-THETA( 1 8 ) 

C PANG(18)-PANG(18)+1 

C TYPE *PANG(18) 

C ENDIF 

C XNR15=XNR12+ ( 3*XNRATA) 

C TYPE *,’XNR15-\XNR15 

C IF (RANDTTJCNR15ANDJRAND.GE.XNR14) THEN 
C PHI3-PHI(19) 

C THETA3 -THETA ( 19 ) 

C PANG( 19)=PANG(19)+1 

C TYPE*,PANG(19) 

C ENDIF 

C XNR16-XNR12+(4*XNRATA) 

C TYPE *,’XNR16«’,XNR16 

C IF (RANDTTJXNR16ANDRAND.GE-XNR15) THEN 
C PHI3-PHI(20) 

CC THETA3-THETA( 20) 

C PANG ( 20 ) -PANG ( 20 ) + 1 

C TYPE * J*ANG (20) 

C ENDIF 

C XNR17-XNR12+(5*XNRATA) 

C TYPE *,’XNR17-’,XNR17 

C IF ( RAND.LT.XNRATAAND.RAND.GE.0.0 ) THEN 
C PHI3-PHI(21) 

C THETA3-THETA(21) 

C PANO(21)-PANG(21)+l 
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C TYPE*J»ANG(21) 

C ENDIF 
C 

q************************************************ 

c **** THREE POINT THERMAL FLARE******* 
£************************************************ 
c 

C IF ( RAND !£.. 28 ) THEN 
C PHI3-PHI(1) 

C THETA3-THETA( 1 ) 

C PAN0(1)-PAN0(1)+1 

C ENDIF 

C IF ( RAND.GT..28 -ANDJRAND.LE..68 ) THEN 
C PHI3»PHI(11) 

C THET A3 -THETA ( 1 1 ) 

C PANG(11)-PANG(11)+1 

C ENDIF 

C IF (RAND.GT..68)THEN 
C PHI3-PHI(21) 

C THETA3 -THETA (21) 

C PANG(21)-PANG(21)+1 

C ENDIF 
C 

C MODEL T AT TIME-0 AT 10 KEV 

C*********************************************************** 

c 

C XYZ-9.98E-4 

C IF ( RAND LE.XYZ) THEN 
C PHI3=PHI(1) 

C THETA3-THETA( 1 ) 

C PANG(1)«PANG(1)+1 

C ENDIF 

C IF (RAND. GT-XYZAND.RAND.LT. ( .999 ) )THEN 
C PHI3-PHI( II) 

C THETA3«THETA( 11) 

C PANG(11)-PANG(11)+1 

C ENDIF 

C IF (RAND.OT.(. 999)) THEN 
C PHI3-PHI(21) 

C THETA3-THETA(21) 

C PANG(21)-PANO(21)+l 

C ENDIF 
C 

******************************************** 
Q******************************************** ****** 

C MODEL NT AT 10 KEV AT 40 SEC 

q******************************************»******* 


XNTOT -532320.0 
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XNRAT1~( 144000./XNTOT) 

XNRAT2= ( 43200. /XNTOT ) +XNRAT1 
XNRAT3- ( 240. /XNTOT ) +XNRAT2 
XNRAT4= ( 960. /XNTOT ) +XNRAT3 
XNRAT5= ( 1920. /XNTOT ) +XNRAT4 
XNRAT6- ( 28 80. /XNTOT ) +XNRAT5 
XNRAT7- (3360. /XNTOT ) +XNRAT6 
XNRAT8- ( 14400. /XNTOT ) +XNRAT7 
XNRAT9- ( 19200. /XNTOT ) +XNRAT8 
XNRAT10- ( 24000. /XNTOT ) +XNRAT9 
XNRAT1 1 -( 24000. /XNTOT ) +XNRAT 10 
XNRAT12- (24000. /XNTOT) +XNRAT1 1 
XNRAT13=( 19200./XNTOT)+XNRAT12 
XNRAT14=( 14400. /XNTOT) +XNRAT1 3 
XNRAT15- (3360./XNTOT) +XNRAT14 
XNRAT16=(2880./XNTOT)+XNRAT15 
XNRAT17- ( 1920. /XNTOT) +XNRAT16 
XNRAT1 8- ( 960. /XNTOT) +XNRAT17 
XNRAT19- ( 240. /XNTOT) +XNRAT1 8 
XNRAT20= (43200. /XNTOT) +XNRAT19 
XNRAT2 1 = ( 144000. /XNTOT ) +XNRAT20 

TYPE *, r RAND-,RAND 
TYPE *, ’XNTOT- ’.XNTOT 
TYPE */XNRATl«\XNRATl 
TYPE *,’XNRAT2-’XNRAT2 
TYPE *,’XNRAT3=’,XNRAT3 
TYPE *,’XNRAT4=’,XNRAT4 
TYPE ♦.’XNRAT5- XNRAT5 
TYPE *.’XNRAT6=\XNRAT6 
TYPE *,’XNRAT7-’^NRAT7 
TYPE *,’XNRAT8=’,XNRAT8 
TYPE VXNRAT9-’,XNRAT9 
TYPE *,’XNRAT10-’^NRAT10 
TYPE *,’XNRAT11«\XNRAT11 
TYPE *,’XNRAT12«\XNRAT12 
TYPE */XNRAT13=’ .XNRAT13 
TYPE *,XNRAT14-’XNRAT14 
TYPE *,’XNRAT15«\XNRAT15 
TYPE VXNRAT16-.XNRAT16 
TYPE *,’XNRAT17-’XNRAT17 
TYPE *,’XNRAT18«\XNRAT18 
TYPE *,*XNRAT19-’,XNRAT19 
TYPE * , ’XNRAT20 - ’ XNRAT20 
TYPE *,’XNRAT21»’,XNRAT21 

IF (RANDiTJCNRATl-ANDHAND.OE.0.0) THEN 


C 

C 


PHI3-PHI(1) 

THETA3-THETA(1) 

PANG(1)«PANG(1)+1 
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C TYPE*,PANG(1) 

C TYPE *,TLARE POSITION- ’,1 

C ENDIF 

C IF ( RAND.LT. ( XNRAT2 ) ANDRANDGEJXNRAT1 ) THEN 
C PHI3-PHI(2) 

C THETA3-THETA(2) 

C PANG(2)-PANG(2)+1 

C TYPE*J»ANG(2) 

C TYPE ♦,’FLARE POSITION -’,2 

C ENDIF 

C IF (RANDTT. (XNRAT3) AND JIAND.GE.(XNRAT2)) THEN 
C PHI3-PHI(3) 

C THETA3-THETA ( 3 ) 

C PANG(3)-PANG(3)+1 

C TYPE*PANG(3) 

C TYPE ♦,’FLARE POSITION-’, 3 

C ENDIF 

C IF (RANDX,T.(XNRAT4) AND JIAND.GE.(XNRAT3)) THEN 
C PHI3-PHI(4) 

C THETA3-THETA ( 4 ) 

C PANG(4)-PANG(4)+1 

C TYPE *J*ANG(4) 

C TYPE ♦,’FLARE POSITION- ’,4 

C ENDIF 
C 

C IF(RAND.LT.(XNRAT5)AND.RAND.GE.(XNRAT4))THEN 
C PHI3-PHK5) 

C THETA3-THETA(5) 

C PANG(5)-PANG(5)+1 

C TYPE*PANG(5) 

C TYPE ♦,’FLARE POSITION-’^ 

C ENDIF 
C 

C IF ( RAND J-T.(XNRAT6) AND AAND.GE.(XNRATS) ) THEN 
C PHI3-PHI(6) 

C THETA3-THETA( 6) 

C PANG(6)-PANG(6)+1 

C TYPE*J»ANG(6) 

C TYPE ♦,’FLARE POSITION-’, 6 

<X ENDIF 
C 

C IF (RANDTT. ( XNRAT7 ) ANDRAND.GE. ( XNRAT6 ) ) THEN 
C PHI3-PHI(7) 

C THETA3-THETA(7) 

C PANG(7)-PANG(7)+1 

C TYPE *JPANG(7) 

C TYPE ♦,’FLARE POSITION- ’,7 

C ENDIF 
C 

C IF (RANDET.(XNRAT8) AND RAND . GE.(XNRAT7) ) THEN 
C PHI3-PHH8) 
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c THETA3=THETA( 8 ) 

C PANG(8)«PANG(8)+1 

C TYPE ‘.PANG ( 8 ) 

C TYPE *,’FLAREPOSmON=\8 

C ENDIF 
C 

C IF (RANDiT.(XNRAT9) AND JIAND.GE.(XNRAT8)) THEN 

C PHI3=PHI(9) 

C THETA3«THETA(9) 

C PANG(9)=PANG(9)+1 

C TYPE *,PANG(9) 

C TYPE ‘/FLARE POSITION-'^ 

C ENDIF 
C 

C IF ( RAND I.TJXNRAT10 AND -RAND.GEXNRAT9) THEN 
C PHI3=PHI( 10) 

C THETA3=THETA( 10) 

C PANG(10)=PANG(10)+1 

C TYPE *,PANG( 10) 

C TYPE ‘.TIARE POSITION-’, 10 

C ENDIF 
C 

C IF (RAND.LT JCNRAT1 1 AND.RAND.GE .XNRAT10) THEN 
C PHI3-PHI(11) 

C THETA3 -THETA ( 1 1 ) 

C PANG(11)-PANG(11)+1 

C TYPE*,PANG(11) 

C TYPE ‘/FLARE POSITION=’,l 1 

C ENDIF 
C 

C IF (RANDJLTJXNRATI2ANDJIAND.GE-XNRAT11) THEN 
C PHI 3 -PHI ( 12) 

C THETA3=THETA( 12) 

C PANG(12)-PANG(12)+I 

C TYPE ‘^ANG (12) 

C TYPE ‘/FLARE POSITION =’.12 

C ENDIF 
C 

C IF (RANDJ-T.XNRAT13 AND.RAND.GE.XNRAT12) THEN 
CC PHI3»PHI( 13) 

C THETA3-THETA( 1 3 ) 

C PANG(13)-PANG(13)+1 

C TYPE*J*ANG(13) 

C TYPE*, ’FLARE POSITION-M3 

C ENDIF 
C 

C IF (RAND.LT.XNRAT14AND.RAND.GE.XNRAT13) THEN 
C PM3=PHI(14) 

C THETA3=THETA( 14) 

C PANG ( 14) -PANG ( 14) + 1 

C TYPE ‘.PANG ( 14) 
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C TYPE *, "FLARE POSITION- ’,14 

C ENDIF 
C 

C IF ( RAND TTJCNRAT15 AND JIAND.GEJCNRAT14) THEN 
C Pffl3-PHI(15) 

C THETA3 -THETA ( 15 ) 

C PANO(15)-PANG(15)+l 

C TYPE*PANG(15) 

C TYPE *, "FLARE POSITION-’, 15 

C ENDIF 
C 

C IF (RANDiTJCNRAT16ANDJlAND.GEJCNRAT15) THEN 
C PHI3-PHI(16) 

C THETA3-IHETA( 16) 

C PANG(16)-PANG(16)+1 

C TYPE *J*ANG( 16) 

C TYPE ♦,’FLARE POSITION-’, 16 

C ENDIF 
C 

C IF (RANDiTJCNRAT17ANDJlAND.GEJCNRAT16) THEN 
C PHI3=PHI(17) 

C THETA3 -THETA ( 17) 

CC PANG(17)-PANG(17)+1 

C TYPE *J*ANO( 17) 

C TYPE*, "FLARE POSITION-", 17 

C ENDIF 
CC 

C IF (RANDTTJCNRAT18ANDJLAND.GE.XNRAT17) THEN 
C PHI3-PHI(18) 

C THETA3-THETA ( 18 ) 

C PANO(18)-PANG(18)+l 

C TYPE *,PANG( 18) 

CC TYPE *, ’FLARE POSITION- ’,18 

CC ENDIF 
C 

C IF (RAND1T JCNRAT19 ANDJIAND.GE JCNRAT18 ) THEN 
C PHI3-PHI(19) 

C THETA3-THETA( 19 ) 

C PANO(19)-PANG(19)+l 

C TYPE*J»ANG(19) 

CC TYPE *, "FLARE POSITION-’, 19 

C ENDIF 
C 

C IF ( RAND.LT JCNRAT20 AND JIAND.GE JCNRAT19 ) THEN 
C PHI3-PHI(20) 

C THETA3 -THETA ( 20 ) 

C PANG (20) -PANG ( 20) +1 

C TYPE *PANQ(20) 

C TYPE *, "FLARE POSITION-’, 20 

C ENDIF 
C 
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C IF ( RANDi-T JCNRAT2 1 AND.RAND.GE.XNRAT20 ) THEN 
C PHI3-PHI(21) 

C THETA3=THETA(21 ) 

C PANG(21)=PANG(21)+1 

C TYPE*FANG(21) 

C TYPE *, ’FLARE POSITION=\21 

C ENDEF 
C 


£***************************************************** 

Qt**************************************************** 


(l************************************************** 


C MODELT AT 10KEV AT 10 SEC 


£************************************************** 


XNRAT1=.0769 

C 

IF (RANDiTJCNRATl ANDJtAND.GE.0.0) THEN 
C 

PHD-PHI(l) 

THETA3-THETA(1) 

PANG(1)*PANG(1)+1 
C TYPE *,PANG( 1 ) 

C TYPE ♦,’FLARE POSITION-’, 1 

ENDIF 
C 

IF (RAND.LT.(2*XNRAT1) AND JIAND.GE JG4RAT1 ) THEN 
PHI3=PHI(6) 

THETA3-THETA(6) 

PANG(6)=PANG(6)+1 
C TYPE *,PANG ( 2 ) 

C TYPE *, ’FLARE POSITION-’^ 

ENDIF 

C 

IF (RAND.LT. ( 3*XNRAT1 ) AND JtAND.GE. ( 2*XNRAT1 ) ) THEN 
PHI3=PHI(7) 

THETA3«THETA(7) 

PANG(7)-PANG(7)+1 
C TYPE*J»ANG(3) 

C TYPE *, FLARE POSITION-’, 3 

ENDIF 
C 

IF (RANDiT.(4*XNRATl)ANDJLAND.GE.(3*XNRATl)) THEN 
PM3«PHI(8) 

THETA3«THETA(8) 

PANG(8)-PANG(8)+1 
C TYPE * J*ANG(4) 

C TYPE VFLARE POSITION-’, 4 

ENDIF 
C 
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IF (RANDIM\(5*XNRATl)ANDJtAND.GE.(4*XNRATl)) THEN 
PHI3-PHI(9) 

THETA3=THETA(9) 

PANO(9)-PANO(9)+l 
C TYPE*PANG(5) 

C TYPE ♦,'FLARE POSITION- ,5 

ENDEF 
C 

IF (RANDiT.(6*XNRATl) ANDJtAND.GE.(5*XNRATl)) THEN 
PHI 3 -PHI ( 10) 

THETA3-THETA( 10) 

PANG(10)-PANG(10)+1 
C TYPE*PANG(6) 

C TYPE *, FLARE POSITION-’, 6 

ENDIF 
C 

IF (RAND.LT. ( 7*XNRAT1 ) AND JIAND.GE. ( 6*XNRAT1 ) ) THEN 
PHI3-PHI(11) 

THETA3-THETA( 11) 

PANG(11)-PANG(11)+1 
C TYPE * J*ANG ( 7 ) 

C TYPE *, FLARE POSITION-’, 7 

ENDIF 
C 

IF ( RAND.LT. ( 8 *XNRAT1 )AND.RAND.GE.( 7*XNRAT1) ) THEN 
PHI3-PHI(12) 

THETA3-THETA( 12) 

PANG(12)-PANG(12)+1 
C TYPE * J*ANG ( 8 ) 

C TYPE ‘.’FLARE POSITION- ’,8 

ENDIF 
C 

IF (RANDiT.(9*XNRATl) ANDJtAND.GE.(8*XNRATl)) THEN 
PHI3«PHI( 13) 

THETA3-THETA( 13) 

PANG(13)-PANG(13)+1 
C TYPE *,PANG(9) 

C TYPE *, FLARE POSITION- ’,9 

ENDIF 
C 

IF (RANDLrT. ( 10*XNRAT1 ) ANDRAND.GE. ( 9*XNRAT1 ) ) THEN 
PHI3-PH3(14) 

THETA3 -THETA ( 14 ) 

PANG ( 14) -PANG ( 14) +1 
C TYPE*J»ANG(10) 

C TYPE*, FLARE POSITION- ’,10 

ENDIF 


C 


IF (RAND.LT.( 1 1*XNRAT1 ) ANDRAND.GE.( 10*XNRAT1 ) ) THEN 
PHI3«PHI(15) 

THETA3-THETA( 15) 
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PANG(15)=PANO(15)+l 
C TYPE*,PANG(11) 

C TYPE *,’FLAREPOSmON=\ll 

ENDIF 
C 

IF (RAND.LT. ( 12*XNRAT1 ) ANDHAND.GE. ( 1 1 *XNRAT1 ) ) THEN 
Pffl3-PHI(16) 

THETA3»THETA( 16) 

PANG(16)-PANO(16)+l 
C TYPE *,PANG( 16) 

C TYPE *, ’FLARE POSITION-M6 

ENDIF 


C 

IF (RANDJLT.( 13*XNRAT1 ) AND JIAND.GE. ( 12*XNRAT1) ) THEN 
PHI3-PHI(21) 

THETA3-THETA(21 ) 

PANG(21)-PANG(21)+1 
C TYPE *J*ANG(21 ) 

C TYPE *,’FLARE POSITION=\21 

ENDIF 


Q***************************************************** 


RETURN 

END 

SUBROUTINE XSTARXY 


CALCULATE IMPACT OF PHOTON ON LOWER GRID 


INCLUDE ’XSTAR.INC/LIST’ 


TYPE *, ’ENTERING XSTARXY’ 

WRITE (3,*) ’ENTERING XSTAR XY’ 

G SEP-520 !CM 

TYPE *,THETA3(SUB)-’,THETA3 
DELX-GSEP*TAN ( PHI3 ) *SIN ( THETA3 ) 
DELY-GSEP*TAN ( PHI3 ) *COS ( THETA3 ) 

XPOSX-XPOS-DELX 

YPOSX-YPOS-DELY 

XPOS2-(XPOSX*COS(TWIST) )+(YPOSX*SIN(TWIST) ) 
YPOS2- ( YPOSX*COS ( TWIST ) ) - ( XPOSX*SIN ( TWIST ) ) 

WRITE( 3 ,*) ’XPOS2-’XPOS2,’YPOS2-’,YPOS2 

IF ( XPOS2.GT. 10.0.OR.YPOS2.GT. 10.0 ) THEN 
PTEST-7 
ENDIF 
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c 

IF (XPOS2LT.O.OOR.YPOS2LT.O.O) THEN 
PTEST-7 
ENDIF 

C TYPE *,’XPOS2-’,XPOS2,’YPOS2-’,YPOS2 
C TYPE ♦.’PTEST-’PTEST 
RETURN 
END 

SUBROUTINE XSTARXY2 

CALCULATE IMPACT OF PHOTON ON LOWER GRID 

INCLUDE *XSTAR JNC /LIST’ 

C 

TYPE * ' ENTERING XSTAR XY2’ 

C 

GSEP2-1.0 !CM 
C 

DELX-GSEP2*TAN ( PHI3 ) *SIN ( THETA3 ) 
DELY-GSEP2*TAN(PHI3)*COS(THETA3) 

C 

XPOS3-XPOS2-DELX 

YPOS3-YPOS2-DELY 

C 

IF ( XPOS3.GT. 10.0.OR YPOS3.GT. 10.0 ) THEN 
PTEST-7 
ENDIF 
C 

IF (XPOS3LT.O.O.OR.YPOS3LT.O.O) THEN 
PTEST-7 
ENDIF 
C 

C TYPE*.THI3-,PHI3,”,THETA3- , ,THETA3 
C TYPE *,T)ELX-’ J)ELX,’DELY- , J)ELY 
TYPE *,’XPOS3«’ ( XPOS3,YPOS3-\YPOS3 
RETURN 
END 
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APPENDIX D 


RMC Code 

This appendix contains an example of the code defining the RMC design 
described in Chapter VI an its associated end-to-end, random photon counting 
simulation. The design/simulation provides a data stream of Fourier components 
which may be used as input for AIPS. This operating system will then provide 
a final image. 
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PROGRAM XRAYS MAIN 


C 

C THIS IS THE MOTHER PROGRAM FOR DETERMINING A TELESCOPE’S RESPONSE 
C TO AN INCOMING EXTENDED SOURCE WAVEFRONT USING THE RAY TRACE 
C METHOD AND WRITING THE VISIBILITY POINTS TO A U,V FILE FOR INPUT 
C INTO AIPS. IN THIS PROGRAM, THE EXTENDED SOURCE IS TREATED AS A 
C COMBINATION OF POINT SOURCES WHICH ARE SUPERIMPOSED TO FORM A 
C COMPLEX TRACE WHICH IS THEN BROKEN INTO COMPONENTS USING THE BINNI 
C METHOD. 

C 

INCLUDE ’FEKESDATAJNC/LIST’ 

TELESCOPE GEOMETRY LOOP 

DO 1000NGRID = U 

TYPE *, ’NGRID - ’JMGRID 

TYPE * ’ ENTERING FIXESGEO’ 

CALL FIXESGEO (NGRID) 

CALL GRID_GEO( NGRID) 

FLARE GEOMETRY LOOP 

DO 1000 K5-1.20 
C 

TYPE ♦,’FLARE GEOMETRY LOOP: ’SOURCE POINT # = ’JK5 
TYPE *, ’ENTERING FLAREGEO’ 

C 

CALL FLARE GEO ( K5 ) 

C 

TYPE *, ’ENTERING XRAYSBETA’ 

C 

CALL XRAYS BETA 
C 

TYPE *, ’ENTERING XRAYSTHARY* 

C 

CALL XRA Y STHAR Y 
C 

TYPE *, ’ENTERING XRAYSTHETA’ 

C 

CALL XRAYS THETA 
C 

TYPE *, ’ENTERING XRAYSSLOPE’ 

C 

CALL XRAYS SLOPE 
C 
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TYPE *, ’ENTERING XRAYSMVAL’ 

C 

CALL XRAYS MVAL 
C 

TYPE *, ’ENTERING XRAYSPHI’ 

C 

CALL XRAYS PHI 
C 

TYPE *, ’ENTERING XRAYSALPHA’ 

C 

CALL XRAYS ALPHA 
C 

TYPE *, ’ENTERING XRAYSF ALPHA’ 

C 

CALL XRAY S_F ALPHA 
C 

TYPE *, ’ENTERING XRAYSFALPHAROTATE’ 

C 

CALL XRAYS FALPHA ROTATE 
C 

TYPE * ’ENTERING XRAYSSAWTH’ 

C 

CALL XRAYS SAWTH ( NGRID,K5 ) 

C 

C900 CONTINUE 
C 

TYPE *, ’ENTERING XRAYSBIN’ 

C 

C NOTE: THIS PLOT CALL MUST BE COMMENTED 
C OUT TO USE THE ONE IN SAWTH 
C 

C IF (K5.EQ. 1 .OR. 20) THEN 

C CALL PLOT XRAY ( SAWTP,FSAWTP > IB,8 ( DVAL,LVAL,KI) 
C ENDIF 
C 

CALL XRAYS_BIN(NGRIDJC5) 

C 

1000 CONTINUE 
C 

9999 STOP 
END 

C FIXESDATA.INC 
C 

REAL*4 THARY(2000) .FTHARY ( 2000) .THETA ( 8100), 

1 FTHETA(8100) J>THETA( 1200),SLOPE( 1200) 

C 

REAL*4 BETA( 100) JDPHH 1000) PHI( 30000) 

C 

REALM SAWT(2, 30000) JFSAWT(2, 30000) 

C 

REALM SAWT1 ( 30000 ) .FSAWTl ( 30000) 


235 


c 

REALM SAWTP ( 30000 ) .FSAWIP ( 30000 ) 

C 

REAL*4 ALPHA ( 30000) F ALPHA ( 30000) ,XPHI( 30000) ,YFALPHA( 30000) 
C 

REALM DVALLVALJ>IMJM( 101 ) L>NUM( 101 ) 

C 

REALM DUM1 ( 30000) ,DUM2( 30000) J3UM3 ( 30000) JXJM4( 30000) 

C 

REALM BINA(200)3INB(200)3INC(200)3IND(200),SUM(200) 

C 

REALM U ( 200 ), V ( 200 ), VX ( 200 ),VY( 200) 

C 

REALM REOA( 200), SAWTH(40),FSAWD(2,30000),TERVAL 
C 

INTEGERM MVAL.NGRID JQ,XL( 8 ) >1 
C 

CCC DATA DVAL/. 0027/LVAL/40.0/ J>I/3. 14159/ JQ/ 10/ 

CCC DATADVAL/.0349/LVAL/40.0/J>l/3.14159/J0/I0/ 

CCC DATA DVAL / .00465 / LVAL/40.0 / J>I/ 3. 14159 / JO / 10/ 

CCC DATA DVAL/. 0087/ LVAL/40.0/ J > I/3.14159/JQ/10/ 

CCC DATADVAL/.0175/LVAL/40.0/J»l/3.14159/J0/I0/ 

C 

DATALVAL/5.2/J»I/3.14159/JO/l/ 

C 

COMMON /THET/ THAR Y JTHAR Y, THETA JTHETAJDTHETA, SLOPE 
C 

COMMON /BEPH/ BETAJM»HIJ»HI 
C 

COMMON /ALPH/ ALPHA J?ALPHA,XPHI,YF ALPHA 
C 

COMMON /RDATA/ DVALLVALJ’INUMJDNUM 
C 

COMMON / IDATA / KIJ4VAL 
C 

COMMON/SAW/ SAWTH.SAWTESAWT.SAWTPESAWTP 
C 

COMMON / SAWW/ SAWT13SAWT1 
C 

COMMON /BIN/ BINA3INB3INC3INDJB,SUM 
C 

COMMON/DUM/ DUM1,DUM2,DUM3,DUM4,VX,VY,U,V 
C 

COMMON /REG / REG A, SAWD 
C 

C 

SUBROUTINE FLAREGEO ( K5 ) 

C 

C THE PURPOSE OF THIS SUBROUTINE IS TO BUILD A FLARE 
C FROM A COLLECTION OF POINT SOURCES, A SIGNIFICANT 
C PROBLEM IS WHAT COLLECTION OF POINTS BEST REPRESENT 
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C A FLARE OF INTEREST. 

C 

C THIS ATTEMPT AT PAINTING A FLARE WILL BE REPRESENTED 
C BYA 100 ARCSECOND LONG SEMICIRCLE WITH 20 POINT 
C SOURCES EQUA- SPACED ALONG THE ARC. THIS MEANS THE 
C RADIUS OF THE FLARE WILL BE 100/PI. 

C 

C THE ORIENTATION OF THE FLARE TO THE VERTICAL AXIS 
C OF THE TELESCOPE WILL BE 45 DEGREES. THE FIRST 
C POINT WILL BE 4 5 DEGREES OF ROTATION FROM THE 
C HORIZONTAL, THE NEXT WILL BE 4.5 + 9.0, THE NEXT 
C 45+9.0+9.0 AND SO FORTH. ESSENTIALLY, DIVIDING 
C 180 DEGRESS BY 20 PIXELS AND PUTTING THE POINT 
C SOURCE IN THE CENTER OF EACH PIXEL. 

C 

INCLUDE TTXESDATA.INC /LIST’ 

C 

DIMENSION PSI(30) 


OPEN(UNIT-3 JTLE-’FLAREOAT’.STATUS-’NEW’ ) 


i 


LET PSI BE THE ANGLE AROUND THE CIRCLE 
LET DELPSI BE THE ITERATION OR 9.0 DEGREES 

PSI(l) = 45*(PI/180.0) 

DO J = 2,20 

PSI( J) =PSI( J- 1 ) + (9.0* (PI/ 180.0) ) 

X TYPE *,’I=’,J,TSI(J)«\PSI(J) 

END DO 

DEFINE ANGLES OF FLARE POINT 

ANGLE FROM AXIS OF ROTATION OF TELESCOPE IN RADIANS 
DON -1,20 

WRITE(3,*) TCI BEFORE ERROR=',KI,’N BEFORE ERROR=’,N 
BETA(KI) -BETA(KI) +( (SIN(PSI(N) ))*( 10.0/ 18.0) 

1 *(1/3600.0)) 


WRITE(3,*) *BETA(KI)- , 3ETA(KI),’K1-\KI 
WRITE(3,*) ’PSI(N) )-\PSI(N),’N«’,N 


ANGLE OF ROTATION OF TELESCOPE ABOUT AXIS OF ROTATION 
TO COALIGN POINT SOURCE AND VERTICAL AXIS OF TELESCOPE 


SAWTH(N)-( (PI/180.0)*45.0)+( (COS(PSI(N) ) ) 
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1 *(10.0/18.0)*(l/3600.0)) 

C WRITE (3,*) ’SAWTH(K5)-\SAWTH,’K5-\K5 
C WRITE(3,*) ’BETA(KI)-’,BETA(KI),’KI-\KI 
ENDDO 
RETURN 
END 

SUBROUTINE XRAYSBETA 
C 

INCLUDE ’FDCESDATAJNC/LIST’ 


THE PURPOSE OF THIS LOOT IS TO SPAN THE TELESCOPE FIELD OF VIEW OF 
APROXIMATELY 10 DEGREES FOR A POINT SOURCE. WE WILL INPUT THE 
ANGLE KI. 

DO 1-1,100 

IF(I£Q. 1) THEN 

BETA(I) = PI/ 1800 

ELSE 

BETA IS ANGLE OF POINT SOURCE FROM AXIS OF TELESCOPE 

BETA(I) -BETA(I-l) + PI/1800 
ENDIF 

ENDDO 

RETURN 

END 

SUBROUTINE XRAYSTHARY 
INCLUDE ’FEKESDATAJNC/LIST* 


1-0 


THIS LOOP IDENTIFIES THE NORMALIZED MAXS AND MINS FOR A POINT 
SOURCE AT AN ANGLE THETA THAT THIS DETECTOR SEES 


C2345678901234567890 

C 

DO N- 1,800 
M-N-l 
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1 = 1 + 1 

c 

C THE EQUATION FOR THETA MAX 
C 

THARY(I) = ATAN( (2 * M * DVAL) / LVAL) 

C 

FTHARY(I) - 1.0 
C 

C DIAGNOSTICS: WRITE THETA MAX VALUES TO ’FOR002DAT, FOR VIEWING. 
C 

CCC WRITE(2,25)I,THARY(I),FTHARY(I) 

C 

1 = 1+1 


c 

C THE EQUATION FOR THETA MIN 
C 

THARY(I) = ATAN( ( (2 * M + 1 ) * DVAL ) / LVAL) 

C 

FTHARY(I) =0.0 
C 

C DIAGNOSTICS: WRITE THETA MIN VALUES TO ’FOR002DAT, FOR VIEWING. 
C 

CCC WRITE(2 > 25)I,THARY(I) ( FTHARY(I) 

20 FORMAT( 1X,F12.6,2X,F12.6) 

25 FORMAT! 1X44,’ )',’THARY=’,F12.6,2X, , FTHARY=’.F 12.6) 

END DO 

TYPE *, ’I=’4,THARY(I) =’,THARY(I) ,’FTHARY (I) =\FTHARY (I) 

9999 RETURN 
END 

SUBROUTINE XRA Y STHET A 
C 

INCLUDE TTXESDATA.INC /LIST’ 


C 



c 

C NOW WE BUILD A SAWTOOTH CURVE USING THE END POINTS CALCULATED 
C BEFORE BUILDING AN FTHETA VS THETA SAWTOOTH 
C 


THETA( 1) =THARY(1) 

FTHETA! 1 )= FTHARY ( 1 ) 

1 = 0 

DON =1,800 !500 

C 

FTHETA(N*10) =FTHARY(N+1) 

THETA!N*10) =THARY(N+1) 

C*** DIAGNOSTICS *** 

C 

C WRITE(8,100)NJTHETA(N),THETA(N) 

100 FORMAT! 1X,I4, , )’,’FTHETA= , ,F15.6,2X,’THETA=\F15.6) 
END DO 
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c 

RETURN 

END 

SUBROUTINE XRAYSSLOPE 
C 

INCLUDE TIXESDATAJNC/LIST’ 

C 

C 

C 

C NUMERATOR AND DENOMINATOR DIAGNOSTIC FOR SLOPE 
C 

NUM(l) - (FTHETA(IO)) - (FTHETA(l)) 

DNUM( 1 ) - (THETA( 10) ) - (THETA( 1 ) ) 

C 

C CALCULATING SLOPE OF EACH LINE SEGMENT 
C 

SLOPE (1 ) - NUM( 1 ) /DNUM( 1 ) 

TYPE *,’SLOPE( 1 ) -’,SLOPE( 1 ) 

n-o 

DO ION- 1,800 
M-N+l 
I - (N+ 1) * 10 
J-N*10 
C 

C FOR 1ST SEGMENT, WE MUST USE 9, FOR 2ND, 10 
C 

CCC TYPE*,’F(I)-\FTHETA(I),’F(J)«’,FTHETA(J), 

CCC * T(I)-’,THETA(I),T(J)-’,THETA(J) 

C 

SLOPE(M) - (FTHETA(I)-FTHETA(J)) 

1 /(THETA(I)-THETA(J)) 

C 

IF(N£Q.1)THEN 

DTHETA(N) - (THETA(N*10)-THETA(N))/9 
C 

ELSE 

DTHETA(N) = (THETA(N*10) -THETA((N-1)*10))/10 
C 

END IF 
C 

DOL-1,9 

n-n+i 

c 

THETA(n+l) -THETA(H)+DTHETA(N) 

C 

IF(N£Q.1)1HEN 

FTHETA(E) -SLOPE(N) *THETA(H) + 1.0 
C 

ELSE 

IF ( SLOPE (N).LT.0.0) FTHETA(II+1) - ( SLOPE(N)*(THETA(II+ 1 ) -THETA 
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1 (10*(N-1)))) +1.0 

IF ( SLOPE (N).OT.O.O) FTHETA(II+1 ) = SLOPE(N) * (IHETA(n+l) - 
1 THETA(10*(N-1))) 

ENDIF 

END DO 


C 

n-(N*io)-i 

c 

C SAVEN VALUE INN1 FOR PLOT PURPOSE 
C 

10 CONTINUE 
N1«(N-1)*10 
NNl-n-(L-l) 

TYPE *,’N1=’ ( N1,’NN1= , ,NN1 
C 

c 

C PLOT THETA VS FTHETA 
C CALL PLOT ROUTINE 
C 

C CALLPLOT_XRAY (THETAJ’THETA^OO.IDVALJLVALJO) 

C 

C DO K - 1,1001 

CCC WRITE(7,100)K,THETA(K),FTHETA(K) 

100 FORMAT( 1XJ4,’) \2XTHETA=’JF10.6^,TTHETA=’,F10.6) 

C END DO 
C 

RETURN 

END 

SUBROUTINE XRAYSMVAL 
C 

INCLUDE TTXESDATA.INC /LIST’ 

C 

REALM RVAL 
C LETS CALCULATE MVAL 
C 

DOMN-U00 

C 

RVAL- (2*MN) * (DVAL/LVAL) * (1.0/(TAN(BETA(KI)))) 

C 

IF(RVAL -OT. 1.0) GOTO 10 
C 

MVAL-MN 

C 

RVAL- ((2*MN) + 1) * (DVAL/LVAL) * (10/(TAN(BETA(KI)))) 
C 

IF(RVAL.GT. 1.0) GOTO 10 
C 

CCC WRTTE(8,100)MNMVALJIVAL 

100 FORMAT ( 1X44,’ ) 1X.MVAL- ^.lX.’RVAL-’.F 1 3.9 ) 
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END DO 


C 

10 CONTINUE 

CCC WRITE( 8,101 )MN,MVAL,RVAL 

101 FORMAT ( / / , IX J4,’ ) 1X,’MVAL=’J4, 1X,’RVAL=\F1 3.9 ) 

TYPE *, ’MN-’MN.’MVAL-’MVAL.’RVAL-’.RVAL 

RETURN 

END 

SUBROUTINE XRAYSPHI 
C 

INCLUDE ’FIXESDATAJNC/LIST’ 

C 

C THIS LOOP REPRESENTS THE TELESCOPE ROTATION 
C DEFINES PHI MAXS & MINS AND DELTA PHI’S 
C 

KK-0 

rr-o 

M-MVAL 
DO ION- 1,550 
C 

C PHI IS THE ANGLE OF ROTATION OF TELESCOPE 
C 

c 

IT-IT+ 1 
C 

C GIVES PHI MAX 
C 

Y -1.0/(TAN(BETA(KI) ) ) 

PHI(IT*10)=ACOS( (2 * M * DVAL) /LVAL * Y) 

FALPHA(IT*10)«1.0 

C 

rr-iT+i 

M-M-l 

IF(M.LT.0)GOTO20 

C 

c 

C GIVES PHI MIN 
C 

Y «1.0/(TAN(BETA(KI) ) ) 

PHI (IT* 10) - ACOS ( ( (2 * M+l )*DVAL) /LVAL * Y) 
FALPHA(IT*10)-0.0 
C 
C 

10 CONTINUE 
20 CONTINUE 
C 

C THIS LOOP IS JUST A DIAGNOSTIC WRITE OF THE FALPHA & PHI ARRAYS 
C 

CCC DO J -1,100 

CCC WRITE(8,100)JJFALPHA(J)J>HI(J)X)PHI(J) 

100 FORMAT ( 1XJ4,’ ) ’.TALPHA-’ JF10.6, 1X,’PHI«\F10.6,1X, 
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1 DPHI-’P12.6) 

CCC END DO 



C 

C BREAK PHI INTO EQUAL SEGMENTS 
C 

DO 30 N= 1,800 
IF(N£Q.1)THEN 
DPHI( 1 ) - (PHI( 10) - PHI( 1 ) ) /9 
C 

ELSE 

DPHI(N) - (PHI(N*10) -PHI( (N-l)*10))/10 
ENDIF 
C 

DOK =1,9 
KK-KK+1 
C 

PHI(KK+1)= PHI(KK) + DPHI(N) 

C 

C 

IF(PHI(KK+1) .GE. 1.570796) GO TO 40 
IF(PHI(KK) .GE. 1370796) GO TO 40 
C 
C 

END DO 
C 

KK= (N*10)-l 

KKPTS - KK - (K-l) !THE # OF POINTS TO PLOT 

30 CONTINUE 
40 CONTINUE 

TYPE *, ’N-’J4-1,KK=’KK,THI(KK)-’JPHI(KK) 

c 

c 

C THIS LOOP IS JUST FOR DIAGNOSTIC PURPOSES..PRINT FALPHA & PHI ARRAYS 
C 

CCC DOM - 1,100 

CCC WRITE(8,110)MJ 7 ALPHA(M)PHI(M) J)PHI(M) 

1 10 FORMAT! 1XJ4,’) ’.’FALPHA*’ J J 10.6,1X,’PHI= , P10.6,1X, 

1 T)PHI- , ,F12.6) 

CCC END DO 

C 

C 

C THIS LOOP IS JUST FOR DIAGNOSTIC PURPOSES . . PRINT FALPHA & PHI ARRAYS 
C 

CCC DOM -2500,3000 

CCC WRITE(8.120)MpALPHA(M)PHI(M) 

120 FORMAT! 1XJ4,')’, , FALPHA='P10.6,1X,THI=’P10.6) 

CCC END DO 

TYPE *,’KKPTS=’KKPTS 
C 

C PLOT PHI VS FALPHA 
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C CALL PLOT ROUTINE 
C 

C CALL PLOTXRAY (PMJ 7 ALPHA > KKPTS,2 > DVAL^VAL > KI) 

C 

RETURN 

END 

SUBROUTINE XRA Y S_ ALPHA 
C 

INCLUDE ’FEXESDATAJNC/LIST’ 

C 

C PICK ONE VALUE OF I AT THIS POINT ,1-10 

C 

C 

DON -1,8000 
C 

C ALPHA IS EQUIVALENT VIEW ANGLE COMBINING PHI AND BETA AND 
C CORRESPONDS TO THETA 
C 
C 
C 

ALPHA(N) =ABS(ATAN(TAN(BETA(KI)) *COS(PHI(N)))) 

C 

CCC WRITE(7,100)NALPHA(N)PHI(N) 

100 FORMAT ( 1X,I4, , )’,1X, , ALPHA-’,F10.6,1X,’PHI=\F10.6) 

C 

IF ( ALPHA(N) LE. 0.0) GO TO 10 
C 

END DO 
10 CONTINUE 
N3-N-1 
C 

TYPEVN3«\N3 

C 

C PLOT ALPHA VS PHI 
C CALL PLOT ROUTINE 
C 

C CALL PLOT XRAY (PHIALPIULN3,3J)VALJ.VALJK3) 

C 

RETURN 

END 

SUBROUTINE XRA Y S_F ALPHA 
C 

INCLUDE TIXESDATAJNC /LIST’ 

C 

REALM EXPKHARY( 8000) 

INTEGER *4 N2 
C 

C THE PURPOSE OF THIS LOOP IS TO IDENTIFY THE SEGMENTS FTHETA VS THETA 
C WHICH ALPHA FALLS AND HENCE CALCULATE FALPHA 
CC WRITE(8.96) 

96 FORMAT(lHl) 
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c 

DO 20N -1,8000 
C 

DO 10 K= 1,800 

C 

C TYPE *, T K«’,K-1,’N=’,N,ALPHA(N)=’ALPHA(N), 

C * THARY(K)=’,THARY(K) 

C 

c 

IF( ALPHA(N) .OT.THARY(K) -AND. 

1 ALPHA(N).LT.THARY(K+1)) THEN I2NDIF 

C 

EXPK-((-1.0)**K) 

IF(EXPK JLT.0.0)THEN !3RD IF 

FALPHA(N)=SLOPE(K)*(ALPHA(N)-THARY(K))+1.0 

C 

PI CP 

EXPK = ((-1.0)**K) 

IF(EXPK.GT.0.0)THEN I4THIF 
FALPHA(N) -SLOPE(K) * ( ALPHA(N) -THARY(K) ) 

ENDIF END OF 4TH IF 

ENDIF END OF 3RD IF 

C 

CCC IF(FALPHA(N) .GT. 0.999) FALPHA(N)=1.0 

CCC IF(FALPHA(N) ET. 0.0006) FALPHA(N)-0.0 

C 

HARY(N) -THARY(K) 

GO TO 20 
C 

ENDIF END OF 2ND IF 

C 
C 

10 CONTINUE 
C 

c 

c 

c 

20 CONTINUE 
C 

C SAVE N VALUE IN N2 FOR PLOT PURPOSE 

C 

N2-N-1 

C 

C DIAGNOSTICS 
C 

CCC DON -1,1200 

CCC WRITE(8,100)NALPHA(N)J } ALPHA(N),SLOPE(N),THARY(N)JIARY(N), 

CCC 1 HARY(N+1) 
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CCC100 FORMAT(1XJ4,’)’,1X,’ALPHA-’J 7 10.6,1X,’FALPHA-’ > F10.6,1X, 
CCC 1 ’SLOPE= ’.F17.6, 1X.THARY = ’ JF10.6, 1X,’HARYN=’,F10.6, IX, 

CCC 2 ’HARYN1»\F10.6) 

C 

CCC END DO 
C 

CCC WRITE(8,101) 

CCC101 FORMAT(lHl) 

C 

CCC DON -1,1200 

CCC WRTIE(8,102)N,THARY(N),HARY(N),HARY(N+1) 

CCC102 FORMAT ( 1XJ4,*)\1X,THARY-\F 10.6,1X,’HARYN-’,F10.6.1X, 
CCC 2 ’HARYN1-\F10.6) 

C 

CCC END DO 
C 

TYPE *,’N2-’N2 
C 

C PLOT PHI VS FALPHA 
C CALL PLOT ROUTINE 
C 

C CALLPLOT_XRAY(PHIfALPHAN2,4DVALi,VALJKI) 

C 

RETURN 

END 

SUBROUTINE XRAY SFALPHAROT ATE 
C 

INCLUDE ’FDCESDATAJNC/LIST’ 

C 

INTEOERMK1 

C 

C THIS PART CALCULATES MIRROR IMAGE OF FIRST 90deg OF ROTATION 
C TAKING ADVANTAGE OF SYMMETRY OF PHYSICAL SITUATION. 

C 

M-0 


DON -1,8000 
C 

M-M+ 1 

Ml - ((2*MVAL) + 1) * 10 

M2 -Ml +M 

M3 - Ml -M 
C 

PHI (M2) - (PHI(Ml) -PHI(M3)) +PHI(M1) 
C 

FALPHA(M2) -FALPHA (M3) 

C 

IF(PHI(M2) .GE. PI) GOTO 10 
C 

END DO 
10 CONTINUE 
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c 

K1-M2 

TYPE *, ’M-’,M,’N-'N, , M2- , .M2, , K1- , ,K1 
C 

C THIS LCK)P IS JUST A DIAGNOSTIC 
C 

CCC WRITE(8,96) 

CCC DO J -1,3000 

CC WRITE(8,100)JPHI(J),FALPHA(J) 

C100 FORMAT ( 1XJ4,’)’,’PHI=’,F15.6,1X,’FALPHA-’,F10.6) 

CCC END DO 
C 

C PLOT PHI VS F ALPHA... THE FIRST 90degs PORTION 
C CALL PLOT ROUTINE 
C 

C CALL PLOTXRAY (Pffl,FALPHA,M2,5,DVALLVAL,KI) 

C 

C PASS ’Kl’ VALUE TO THE FALPHA ROTATION SUBROUTINE. K1 IS THE 
C LAST POINT LOCATION IN THE FIRST 90degs PORTION. 

C 

CALL XRAY SPLOTFALPHI ( K 1 ) 

C 

RETURN 

END 

SUBROUTINE XRAYSSAWTH (NGRID.K5 ) 

C 

C THIS SUBROUTINE ROTATES THE COMPRESSED SAWTOOTH TRACES 

C THRU THE APPROPRIATE ANGLE AND THEN ADDS THEM TOGETHER 
C ONE SHOULD HAVE SIX SINUSOID CURVES AT THE END, ONE FOR 

C EACH GRID, THEN THESE ARE BINNED TO GET FOURIER COMPONENTS 
C 

INCLUDE ’FIXESDATA.INC /LIST’ 

C 

C LOOP TO SHIFT SAWTOOTH BY ANGLE SAWTH 
C 

M-NGRID 

C 

C TYPE*,’ NOTE: INCLUDE FILE: (MN)“( 6,24000)’ 

C TYPE *,’ CORRESPONDS TO NGRID - 6’ 

C 

CLOOP 

C 

DON- 1,29080 

SAWT1 (N) =PHI(N)+SAWTH(K5) 

FSAWT1 (N)-FALPHA(N) 

END DO 
C 

C SOME DIAGNOSTICS HERE 
C 

CCC DO NO- 1,1000 

CCC WRITE( 7,* ) NO, ’ SAWTl(NO)=’,SAWTl(NO),'FSAWTl(NO)=’J’SAWTl(NO) 
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CCC * ,’FALPHA(NO ) = ’J 7 ALPHA ( NO ) ,’SAWTH (K5 ) = \SAWTH(K5 ) 
CCC END DO 
C 
C 

C 2D LOOP 
C 

DO J= 1,29080 

IF (SAWT1(J).GE.(2*PI))THEN 
JA-J 
GO TO 10 
ENDIF 
END DO 
C 

10 CONTINUE 

C WRITE(7,*) ’JA-’JA,’ SAWT1(J)-’,SAWT1(I) 

C WRITE(7*) , J- , J 
C 

C 3D LOOP 
C 

DO K« 1,29080 
KB-JA+K-1 

SAWT(MJK)-SAWT1(KB)-(2*PI) 

C 

C 4/4/90 CHANGE HERE 
C 

C FSAWT ( MX ) -FSAWT(MJC ) +FSAWT1 (KB) 

C 

FSAWT (MX ) -FSAWT 1 (KB) 

C 

IF(SAWT(MX)-GE.SAWTH(K5) )1HEN 
IA-K 
GO TO 20 
ENDIF 
END DO 
C 

20 CONTINUE 

C WRITE (7,*) ’IA-’JA,’ KB-’JCB 
C WRITE(7 ,*) ’K-’X,’SAWT(MX)-’.SAWT(MX) 

C WRITE(7 *) ’M-’JM,’ SAWTH(K5)-’,SAWTH(K5) 

C 

C 4TH LOOP 
C 

DO 1-1,29080 
IC-I+IA 

SAWT(M4C) -SAWT1 ( I ) 

C 

C 4/4/90 CHANGE HERE 
C 

C FSAWT ( MJC ) -FSAWT ( M,IC ) +FSAWT1 ( I ) 

C 

FSAWT(MJC)-FSAWT1(I) 
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c 

IF(SAWT(MJC).GE.(2*PI))OOTO 30 
END DO 
C 

30 CONTINUE 

C WRITE(7 ,*) ’IC-’JC,’ SAWT(MJC)-\SAWT(M,IC) 
C WRITE(7,*) ’FSAWT(M,IC)«’,FSAWT(M,IC) 

C 

C5THLOOP 

C 

DO IB-1,29080 
SAWTP(IB)-SAWT(MJB) 

FSAWTP(IB) -FSAWT(MJB) 

IF (SAWTP(IB).GE.(2*PI))GOTO40 
END DO 
40 CONTINUE 


IG-IB 

DOIF-1,10 

IG-IG-1 

WRITE (7,*) ’NGRID-’JMGRID,’ TM-’M.’ ’.’N-’.N 
WRITE(7 ,*) TB-’JB,’ IG-MG,’ SAWTH(K5)-\SAWTH(K5) 
WRITE(7,*) ’IG-1=’,IG-1,’SAWTP(IG-1)=\SAWTP(IG-1) 
WRITE(7,*) , SAWTP(IG-2)-’,SAWTP(IG-2) 

WRITE(7*) TSAWTP(IG-1)-\FSAWTP(IG-1) 

WRITE(7,*) FSAWTP(IG-2)-’,FSAWTP(IG-2) 

WRITE(7,*) ’SAWT (MJG- 1 ) =\SAWT (MJG- 1 ) 

WRITE(7 ,*) , SAWT(MJG-2)-’,SAWT (M.IG-2) 

WRITE(7 ,*) ’FSAWT(M,IG-1)-\FSAWT(M,IG-1) 

WRITE(7 ,*) TSAWT (M.IG-2 ) - ’ JFSAWT ( MJG -2 ) 

END DO 

IF (K5£Q.20)THEN 

CALL PLOT XRAY ( SAWTPFSAWTP,IB,8 > DVALLVAL,KI) 
ENDIF 

DIAGNOSTIC WRITE 

RETURN 

END 

SUBROUTINE XRAYSBIN (NGRTOJC5 ) 

THE PURPOSE OF THIS SUBROUTINE IS TO BIN THE 
FINAL SINUSOIDAL FUNCTION FOR EACH GRID AND 
PRINT IT OUT TO AN FORTRAN OUTPUT FILE 

INCLUDE ’FDOESDATAJNC/LIST’ 
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c 

C INITIALIZATION 
C 

CHK-0.0 
TERVAL4 - 0.0 
J-I 

CCCC SUM-0.0 
1-0 
C 

M-NGRID 

C 

C TYPE * .’NOTICE HERE THAT XL(NGRID) MAY NEED TO’ 
C TYPE *,TJE VARIED TO OPTIMIZE TELESCOPE’ 
c 

XL(1)=32 !8 

XL(2)-84 !16 

C XL(3)-32 
C XL(4)-64 
C XL(5)«128 
C XL(6)-128 
C 

TERVAL - PI/ (4*XL(M) ) 

C 

DUMKD-0.0 

DUM2(l)-0.0 

DUM3(1)«0.0 

DUM4(l)-0.0 

C 

C IF (K5.EQ.1 ANDNGRIDLQ. 1 ) THEN 

OPEN ( UNIT-4 .FILE- ’IN 2FEL .DAT’, STATUS -NEW’) 

C ENDIF 
C 

c 

c 

c 

C LOOP FOR 360 DEGREES 
C 

DO WHILE (SAWT(MJ)LE.(2*PI)) 

C 

C INITIALIZATION 
C 

IF (K5£Q.l) THEN 
BINA(I)-0.0 
BINB(I)-0.0 
BINC(I)-0.0 
BIND (I) -0.0 
ENDIF 
C 

I-I+l 

C 

TERVAL 1 - TERVAL + TERVAL4 
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TERVAL2 - TERVAL + TERVAL1 
TERVAL3 - TERVAL + TERVAL2 
TERVAL4 - TERVAL + TERVAL3 
C 

REGA(I) =TERVAL2 
C 

CCC WRTrE(4,*)’REGA(I)= > ,RE0A(I), ’I-M 
CCC WRITE(4 *) , SAWT(MJ)-’,SAWT(MJ),’ 
TYPE *,’REGA(I)«\REGA(I),’I=\I 
C 

C TYPE 1ST BIN LOOP’ 

C 

C LOOP FOR 1ST BIN 
C 

DO WHILE (SAWT(MJ)LE.TERVALl) 
CCCC IF (FALPHA(J)-EQ. (0.0)) THEN 
BINA(I)=BINA(I)+1.0 
CCCC ENDIF 
J=J+1 
END DO 
C 

C TYPE *,’2D BIN LOOP’ 

C 

C LOOP FOR 2D BIN 
C 

c 

DO WHILE ( SAWT(M J ) LE.TERVAL2 ) 
CCCC IF (FALPHA(J)EQ.(O.O) ) THEN 
BINB(I)«BINB(I)+1.0 
CCCC ENDIF 
J=J+1 
END DO 
C 

C USE MIDPOINT FOR U,V POINT 
C 

C TYPE *,’3D BIN LOOP’ 

C 

C LOOP FOR 3D BIN 
C 

DO WHILE ( SAWT(MJ ) LE.TERVAL3 ) 
CCCC IF (FALPHA(J)EQ.(O.O)) THEN 
BINC(I)«BINC(I)+1.0 
CCCC ENDIF 
J-J+l 
END DO 
C 

C TYPE *,’4TH BIN LOOP’ 

C 

C 

C LOOP FOR 4TH BIN 
C 
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DO WHILE ( SAWT (MJ) LE.TERVAL4 ) 

CCCC IF (FALPHA(J)3Q.(0.0)) THEN 

BIND(I)-BIND(I)+1.0 
CCCC ENDIF 

J-J+l 
END DO 
C 
C 

TYPE *,’J«\J 

CCC WRITE(4*)’M- , > M, , J-’J 

CCC WRITE (4* ) ’SAWT-\SAWT(M J ) ,’FSAWT(M,J)=’3SAWT(M,J ) 

C 

C CALCULATE U,V 
C 

CCC WRTIE(4,*)’BINA(I)-’3INA(I),’ M-’J 

CCC WRITE(4,*)’BINB(I)- , 3INB(I),’M- , J 

CCC WRTIE(4*)BINC(I)-’3INC(I),’7I-\I 

CCC WRITE(4*)BIND(I)- , 3IND(I), , VI- , ( I 

C 

VX(I)=BINA(I)-BINB(I)-BINC(I)+BIND(I) 

VY(I)-BINA(I)+BINB(I)-BINC(I)-BIND(I) 

U(I)-COS(REOA(I) ) / (DVAL/LVAL) 

V(I)»SIN(REOA(I) ) / (DVAL/LVAL) 

C 

WT-1.0 

C 

C POSSIBLE DIAGNOSTICS 
C 

CHK-CHK+BINA(I)-BINB(I)+BINC(I)-BIND(I) 

SUM(I)«BINA(I)+BINB(I)+BINC(I)+BIND(I) 

JC 

TYPE *,U(I),V(I),VX(I),VY(I) 

C 

C VX(I)»VX(I)*2.0E1 

C VY(I)-VY(I)*2.0E1 

C 

C U(I)-U(I)/25E05 

U(I)— U(I)/2.5E05 
C U(I)-U(I)/2.5E04 

C 

V(I)-V(I)/25E05 
C V(I)-V(I)/25E04 

C 

TYPE *,U(I),V(I), VX(I),VY(I) 

C 

IF (K5EQ.20) THEN 

WRITE ( 4 300 1 ) U ( I ) , V ( I ) , VX ( I ) , VY ( I ) , WT 
2001 FORMAT (2X310.6,1X310.63X311.6,1X310.6,1X310.4) 

ENDIF 

C 

END DO 
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c 

IP-I-1 

c 

IF (K5 EQ. 20) THEN 

CALL PLOTXRAY (REGA,SUMJP,9L)VALLVAL ) KI) 
ENDIF 
C 

TYPE ♦,’NO. OF BINS -’4,’ ,’CHK-’,CHK 
TYPE *, ’SUM- ’.SUM 
C 

RETURN 

END 

X PROGRAM XROT 


*************************************************************** 

** THIS PROGRAM WILL CREATE THE DATA FILE FOR THE ROTATING ** 
****************** MODULATION COLLIMATOR ******************** 
*************************************************************** 


DIMENSION ANGLE(4000) ,U(4000) ,V(4000)X(4000)4a(4000) 
DIMENSION IPASS (4 ),xmin( 25 ),xmax( 25) 

C ** D-Distance between grids ** 

D-500. 

PI-3.14159 

ISEED-3703 

C ** SWID-slit width ** 

SWID1-.0125 

SWID2-.0275 

C * * DIV-number of bins * * 

DIV-1000. 

C ** XMAXPI= maximum value of theta ** 

XMAXPI=2.*PI 


TYPE *, ’Enter angle along x-axis of flare center, in arcseconds’ 
READ *XDELTA 
XDELTA-XDELTA/3600./ 180 *PI 
TYPE *, 'Enter radius of flare, in arcseconds’ 

READ *XPSI 

XPSI-XPSI / 3600. / 1 80.*PI 

TYPE * .Enter number of initial photons’ 

READMTIMES 

TYPE * .Enter number of point sources’ 

read*xflare 

TYPE Enter twist angle, in arcseconds’ 

READ * .TWIST 
TWIST-TWIST/ 3600./ 180.*PI 
type * .Enter number of noise photons’ 
read * .noise 


^ *************************************************************** 
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C ********************** INITIALIZE STUFF ********************* 
q *************************************************************** 


C ** IGRID is the number of the grid set** 

C ** IGRID- 1 is grid set with slit width- 0125 and not shifted ** 

C ** IGRID-2 is grid set with slit width-,0125 and shifted pi/2 * 

C ** IGRID-3 is grid set with slit width=.0275 and not shifted ** 

C ** IGRID-4 is grid set with slit width-.0275 and shifted pi/2 * 

DOIGRID-1,4 

IF (IGR1DBQ.1 .OR. IGRIDJEQ.2) SWID-SWID1 
IF (IGRID.EQ.3 .OR.IGRIDJBQ.4) SWID-SWID2 
IF (IGRIDBQ.l .OR. IGRIDEQ.3) PHASE-0 
IF (IGRID.EQ.2. OR. IGRIDJEQ.4) PHASE-SWID*1.5 
DO I- 1 +DIV* ( IGRID- 1 ) J)IV*IGRID 
ANGLE(I)-0 
U(I)-0 
V(I)-0 
R(I)-0 
XI(I)-0 
ENDDO 

C *************************************************************** 

c ************************** main loop ************************ 

c *************************************************************** 

DO I- 1 +ITIMES* ( IGRID- 1 ) JTIMES*IGRID 

c ******************** FINDS COORD. ON GRID 1 ***************** 

X1-RAN(ISEED)*10. 

Y1-RAN(ISEED)*10. 

Q ***************** CHECKS IF PASSES THRU GRID 1 ************** 
PASS-0 

IDUM1-INT(X1 /SWID) 

IF (IDUMl/2.XEJNT(IDUMl/2.) ) PASS-1 

c ***************** TRANSFORMATION ONTO GRID 2 **************** 

C ** del theta is the angle rotation increment of the collimator * 

IF (PASS.EQ.1 ) THEN 
DELTHETA-XMAXPI /DIV 
THETA-RAN ( ISEED ) *XMAXPI 
C ** dumm is number between 0 & number of point sources ** 
dumm-INT(ran( iseed) *XFLARE ) 

C ** zeta is # between 0 & pi in increments - # point sources ** 

IF (XFLARE£Q.l) ZETA-0 

IF ( XFLARENE. 1 ) ZETA-DUMM*PI / ( XFLARE- 1 . ) 

C ** delta is angle along x-axis of one of the point sources ** 


254 


C ** pgi is angle along y— axis of one of the point sources ** 
DELTA-XDELTA+XPSI*sin(ZETA) 

PSI-XPSI*cos(ZETA) 

C ** sigma the angle of rotation + displacement of point source ** 

C **** due to rotation** 

IF ( DELTA£Q.O ) SIGMA-THETA+PI / 2. 

IF (DELTAHE.O) SIOMA=THETA+ATAN(TAN(PSI) /TAN (DELTA)) 

C * * alpha is displacement of point source along x-axis * * 

C ** beta is displacement of point source along y-axis** 

ALPHA=TAN( SQRT(PSI**2+DELTA**2) ) *COS(SIGMA) 
BETA=TAN(SQRT(PSI**2+DELTA**2) )*SIN(SIGMA) 

X2-X1-D* (ALPHA) 

Y2-Y1 -D*BETA 

q ************************** TWIST *************************** 

C * * finds new values of x2 and y2 on grid 2 due to twist * * 
aa-atan(y2/x2) 
n-sqrt ( x2* *2+y2* *2 ) 
ddd= ( tan ( aa-twist ) ) * *2 
X2=rr/ (sqrt(ddd*( 1+1 / ddd) ) ) 
y2=sqrt(abs(rr**2-x2**2) ) 

q ***************** CHECKS IF PASSES THRU GRID 2 ************** 


X2=X2+PHASE 

IF (X2LE.10 AND. Y2LE.10.) THEN 
IF (X2.GE.0 AND. Y2.GE.0) THEN 
IDUM2-INT(X2/SWID) 

IF (IDUM2/2.NEJNT(IDUM2/2.) ) THEN 
DUM-INT(THETA/DELTHETA) 

IF (DUMNE.O) THEN 
IPASS(IGRID)=IPASS(IGRID)+1 


c **************** ANGLE SUMS PHOTONS THAT PASS *************** 


DUMMY-DUM+DIV* (IGRID- 1 ) 

ANGLE (DUMMY) -ANGLE (DUMMY) + 1 
ENDIF 
ENDIF 
ENDIF 
ENDIF 
ENDIF 
ENDDO 


^ ************************** noise **************************** 


do in- 1 .noise /4 

ntheta-int ( ran( iseed )* 1000 )+ 1 +( igrid- 1 )* 1000 
type *, ntheta 

angle(ntheta) -angle ( ntheta) +1 
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enddo 


c *************************************************************** 
c *************** PREPARES THINGS FOR OUIPUTFILE ************** 
c *************************************************************** 

C ** CORRECTS FOR MISSING AREA OF DETECTOR AT VARIOUS ANGLES ** 

doi-l,div 
theta- (2*pi/div)*i 
dx«abs(D*tan(delta*sin( theta) ) ) 
dy-abs(D*tan(delta*cos(theta) ) ) 
perarea- ( 10-dx) * ( 10-dy ) /100. 
angle(i+ (igrid-1 )*div)=angle(i+(igrid-l )*div) /perarea 
enddo 

C *************** AVERAGES 20 MAX’S AND 20 MIN’S ************** 

do i- 1,25 
xmin(i)« 100000. 
xmax(i)-0. 
enddo 

DO I- 1 +DIV* ( IGRID- 1 ) J>IV*IGRID 
iave-0 
true-1 

dowhile (true.eq.l ) 
iave-iave+1 
if (iaye.gt.24) trne-0 
IF (ANGLE(I).GTXMAX(iave)) THEN 
XMAX(iave)-ANOLE(I) 
true-0 
endif 
en d d o 
iave-0 
true2-l 

dowhile ( true2.eq.l ) 
iave-iave+1 
if (iave.gt.24) true2-0 
if (angle(I) Jtjmin(iave)) then 
xmin ( iave ) -angle ( i ) 
true2-0 

wiHif 

enddo 

enddo 

rain -Q 

smax-0 

doj-2,24 

smin«smin+ rain (j) 
smax~smax+xmax ( j ) 
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enddo 


xave=(smin/24+smax/24) /2 

q ***************** FD'JDSU.V, and real values *************** 


TH-0 

DO I- 1 +DIV* ( IGRID- 1 ) DIV*IGRID 
TH-TH+DELTHETA 
U(I)— D/(2.*SWID)*COS(TH)/10000. 

V(I) =D/ ( 2 *SWID) *SIN(TH) / 10000. 

R(I)-ANGLE(I)-XAVE 

XI(I)=R(I) 

ENDDO 

ENDDO 

q **************** WRITES FILE IN AIPS FORMAT ***************** 


OPEN (UNIT= lJTLE=’XROT DAT’, STATUS=’NEW’ ) 

DO 1-1,1000 

ANGLE(I)=ABS(R(I) )+ABS(XI(I+1000) ) 

WRITE( 1,25) U(I),V(I)^(I)^I(I+ 1000) ANGLE(I) 
ENDDO 
DO 1-2001,3000 

ANGLE(I)-ABS(R(I))+ABS(XI(I+1000)) 

WRTIE(U5)U(I),V(I)Jl(I)^I(I+1000)ANGLE(I) 

ENDDO 


q **************** WRITES FILE WITH PARAMETERS **************** 


OPEN (UNIT-2 J 7 ILE=’XROT2DAT’, STATUS-NEW ) 

WRITE (2,* )DIV,XMAXPI,ITIMES,XDELTA,XPSI,IPASS( 1 ) JPASS(2), 

& IP ASS ( 3 ) JPASS (4) XFLARE,twist,noise 

C ******************** WRITES IDL FILE ************************ 


OPEN (UNIT-3 JFILE=’XROT3DAT, STATUS-NEW’ ) 
WRTTE(3,*)r 

25 FORMAT(2XE10.6,1XE10.6,2XE12.6,1XE12.6,1XE10.4) 

CLOSE(l) 

CLOSE(2) 

CLOSE ( 3 ) 

CALL EXIT 
END 
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